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IV 


CHARACTERIZATION  AND 
FABRICATION  OF  SYNTHETIC  ROUGH 
SURFACES  FOR  ACOUSTICAL 
SCALE-MODEL  EXPERIMENTS 


I.  INTRODUCTION 

The  study  of  energy  scattered  from  rough  surfaces  spans  a  number  of  disciplines,  among  them  those 
involving  either  sound  or  electromagnetic  waves.  Within  the  field  of  acoustics,  rough-surface  scattering 
impacts  sound  fields  in  rooms,  noise  in  urban  environments,  and  acoustic  boundary  interaction  in  the 
ocean.  Ocean  acoustics,  in  particular,  utilizes  theories  of  scattering  from  rough  surfaces  to  characterize  the 
interaction  of  sound  with  the  ocean  surface  (air-water  interface)  and  bottom  (water-bottom  interface  and 
subsurface  layers).  For  characterizing  the  performance  of  sonar  systems,  a  number  of  theories  have  been 
applied  to  cover  a  wide  range  of  frequencies  and  the  entire  range  of  scattering  geometries  in  three  dimensions. 
Some  theories  require  limitations  on  the  roughness,  such  as  perturbation  theory  or  small-slope  theory  [1,  2]; 
others  place  limitations  on  the  geoacoustic  parameters  of  an  ocean  bottom,  such  as  wedge  assemblage  [3,  4]. 
While  rough  surfaces  of  ocean  bottoms  may  involve  topographical  variations  arising  from  either  geological 
processes  or  discrete  objects  (such  as  shells),  here  discussion  focuses  on  rough  surfaces  of  purely  geological 
origin. 

This  report  describes  a  portion  of  a  larger  effort  to  characterize  wave-field  scattering  by  rough  surfaces 
found  in  nature — specifically  the  scattering  of  acoustic  fields  by  rough  ocean  bottoms. 

Characterization  of  acoustic  scattering  by  rough  ocean  bottoms  requires  not  only  devising  an  adequate 
mathematical/physical  theory  of  rough-surface  scattering  and  implementing  it  in  software  (see  [2]),  but 
also  conducting  an  exhaustive  verification  by  comparison  with  scattering  data  from  experimental  measure¬ 
ments  or  numerical  simulations.  However,  both  of  these  methods  for  verification  of  model  predictions  have 
practical  and  technical  limitations.  The  scarcity  of  adequate  data  on  ocean-bottom  properties  and  the  con¬ 
sequent  inability  to  conduct  systematic  studies  of  the  dependence  of  ocean-bottom  scattering  strength  on 
a  roughness  parameter,  such  as  spectral  exponent  (see  Sec.  II B)  make  verification  using  scattering  data 
taken  from  natural  sea  floors  impractical.  Numerical  modeling  (e.g.,  finite-element  and  boundary-element 
methods)  also  has  practical  limitations.  While  it  has  been  successfully  applied  to  two-dimensional  scatter¬ 
ing  from  surfaces  with  stochastic  roughness  that  varies  in  one  dimension  only,  three-dimensional  scattering 
from  surfaces  with  stochastic  roughness  that  varies  in  two  dimensions  can  present  significant  difficulties.  To 
circumvent  these  limitations,  small-scale  physical  models  of  rough  ocean  bottoms  have  been  manufactured 
at  sizes  manageable  for  use  in  an  acoustic  tank  facility  and  verification  has  been  conducted  using  data  taken 
in  the  laboratory  at  appropriately  scaled-up  frequencies.  Such  models  can  be  used  to  investigate  theoretical 
predictions  relating  to  both  roughness  parameters  and  geoacoustic  parameters  such  as  compressional  and 
shear  speeds/attenuations  (see  Ref.  5  and  references  therein).  In  manufacturing  these  small-scale  experi¬ 
mental  models,  no  attempt  has  been  made  to  recreate  any  actual  geophysical  region  in  miniature.  Rather, 
ocean-bottom  topography  is  described  using  the  mathematics  of  self-affine  stochastic  fractals. 

Scattering  theories  applied  in  underwater  acoustics  fall  into  two  categories:  those  involving  a  stochastic 
description  of  the  rough  surface  and  those  treating  the  rough  surface  deterministically.  In  this  report, 
stochastic  rough  surfaces  are  discussed  exclusively.  However,  the  stochastic  rough  surface  is  an  abstraction; 
statistical  properties  of  rough  surfaces  are  more  properly  regarded  not  as  properties  of  individual  surfaces 
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but  as  properties  of  the  ensemble  of  rough  surfaces  that  are  realizations  of  a  particular  random  process. 
Individual  realizations  are,  when  considered  in  isolation,  deterministic.  Therefore,  physical  models  created 
to  represent  stochastic  ocean  bottoms  and  test  stochastic  scattering  theories  can  also  be  used  for  testing 
deterministic  scattering  theories. 

Choosing  to  verify  scattering  theories  by  means  of  physical  scale  models  based  on  fractal  geometry 
requires  that  three  tasks  be  achieved:  (1)  understanding  the  application  of  the  abstract  mathematics  of 
fractal  geometry  to  the  description  of  surface  topography,  (2)  incorporating  into  fractal  models  recent  results 
that  quantify  the  effects  of  representing  unbounded  self- affine  fractal  surfaces  as  bounded,  discretely  sampled 
objects,  and  (3)  developing  or  adapting  the  algorithms  necessary  to  fabricate  physical  realizations  of  fractals 
using  a  modern  computer-numerically  controlled  milling  machine.  The  information  that  was  found  to  be 
essential  for  success  in  these  three  tasks  has  been  summarized  in  this  report.  The  implementation  of 
those  concepts  and  techniques  needed  for  the  manufacture  of  random  rough  surfaces  are  illustrated  by 
accompanying  matlab  code.  This  code  allows  users  to  numerically  generate  rough  surfaces  that  are  two- 
dimensional  realizations  of  band-limited  fractional  Brownian  motion  and  allows  for  such  surfaces  to  be 
characterized  in  terms  of  statistical  and  fractal  metrics. 

This  report  is  grouped  into  two  primary  sections:  Characterization  of  Rough  Surfaces  and  Fabrication  of 
Rough  Surfaces.  In  the  first  of  these  sections,  a  number  of  disparate  elements  are  brought  together  in  order 
to  provide  a  mathematical  description  of  rough  surfaces  as  both  random  processes  and  stochastic  fractals. 
A  parameter  study  (portions  of  which  are  shown  in  the  Appendix  A)  illustrates  the  relationship  between 
the  appearance  of  rough  surfaces  and  their  roughness  parameters.  The  final  portion  of  the  first  section 
relates  the  preceding  theoretical  development  to  real  ocean  bottoms  by  considering  both  how  roughness 
parameters  affect  scattering  and  the  range  of  roughness  parameters  observed  in  the  ocean.  The  second 
section  describes  how  physical  models  of  the  surfaces  described  in  the  first  section  can  be  computationally 
generated  and  then  physically  fabricated  with  the  use  of  computer-aided-manufacturing  equipment,  and  how 
the  accuracy  of  the  resulting  realization  can  be  tested.  The  report  concludes  with  a  summary  of  the  material 
presented.  Appendices  B  and  C  give  further  details  of  the  algorithms  used  for  the  computational  generation 
and  physical  fabrication  of  rough  surfaces,  including  descriptions  of  matlab  software  tools  included  on  the 
CD-ROM  with  this  report.  Abbreviations  and  symbols  are  defined  in  a  glossary  located  at  the  end  of  this 
report.  Numerical  quantities  are  given  in  SI  (MKS)  units  exclusively. 

While  the  primary  goal  of  this  work  is  to  bring  together  concepts  from  disparate  disciplines,  throughout 
the  report  a  number  of  new  concepts  are  introduced.  In  particular,  the  recent  concept  of  approximately 
self-affine  surfaces  is  applied  to  numerical  generation  and  physical  fabrication  of  realizations  of  random  rough 
surfaces.  The  high-  and  low -k  bandwidth  restrictions  inherently  imposed  by  these  processes  are  related  to 
their  effects  on  the  properties  of  surfaces  both  as  characterized  by  fractal  mathematics  and  observed  visually. 


II.  CHARACTERIZATION  OF  ROUGH  SURFACES 

In  this  section,  a  consistent  description  of  rough-surface  topography  is  developed  through  a  series  of 
definitions.  Rough  surfaces  are  described  first  as  infinite-extent  realizations  of  random  processes.  In  Sec. 
II B  each  such  process  is  characterized  by  a  power-spectral  density  (PSD),  in  particular,  a  PSD  that  follows 
a  power-law  form  over  some  portion  of  the  spectrum.  A  parallel  description  of  rough-surface  topography 
in  terms  of  stochastic  fractals  is  developed  in  Sec.  II C  and  related  to  characterization  in  terms  of  PSD. 
Particular  attention  in  both  descriptions  of  surface  roughness  is  given  to  surfaces  having  cutoff  power-law 
PSD  for  which  the  power-law  behavior  is  restricted  in  bandwidth  at  one  or  both  frequency  extremum.  In 
Sec.  II D,  realizations  of  stochastic  rough  surfaces  are  generated  to  provide  an  intuitive  presentation  of  the 
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effects  of  surface-roughness  parameters  on  the  appearance  of  rough  surfaces. 


A.  Statistics  of  Rough-Surface  Topography 

A  rough  surface  is  described  by  its  topography  z  (r) ,  which  is  the  height  of  the  surface  above  a  plane 
containing  position  vector  r  having  Cartesian  components  (x,y).  The  topography  is  given  by  the  mapping 
{z  (r)  :  r  £  R2}  which  defines  a  surface  in  R3  (it  is  assumed  here  that  the  surface  topography  is  single 
valued).  In  this  initial  discussion,  rough  surfaces  are  considered  to  be  of  infinite  extent  and  defined  for  all 
values  of  r,  (i.e.,  ‘continuous’  in  signal-processing  terminology),  as  the  expression  for  the  mapping  from  R2 
implies.  In  practice,  z  (r)  is  defined  over  a  finite  set  of  discrete  points  having  finite  physical  dimensions.  The 
consequences  of  this  will  be  addressed  in  Sec.  Ill  A  1. 

Geomorphological  rough  surfaces  are,  by  necessity  [6],  modeled  as  realizations  of  random  processes  such 
that  the  ensemble  of  realizations  {z}  of  a  particular  random  process  is  described  by  the  probability  density 
function  (PDF)  of  surface  height  fz(()  where  (  is  the  set  of  possible  surface  heights.  In  the  following 
discussion,  it  is  initially  assumed  that  the  random  process  possesses  certain  statistical  attributes  (see,  e.g., 
Ref.  7,  pp.  23-27):  homogeneity  (stationarity  or  translational  invariance),  isotropy  (rotational  invariance), 
and  ergodicity  (equivalence  of  spatial  averages  over  a  single  realization  to  ensemble  averages). 

Given  an  appropriate  choice  of  the  reference  plane  from  which  z  is  measured,  the  mean  height  of  the 
surface  (z)  is  zero 

/oo 

Cfz  (CR  =  o ,  (l) 

-oo 

where  (•)  is  an  expectation  operator.  The  assumption  of  ergodicity  allows  this  operator  to  be  interpreted 
either  as  a  spatial,  (*)s  ,  or  ensemble,  (*)e  ,  average.  Zero- mean  height  can  be  alternatively  achieved  by 
modeling  the  topography  as  the  sum  of  two  functions  z  (r)  +  zq  (r)  where  z  (r)  is  the  stochastic  small-scale 
roughness  and  zq  (r)  is  the  deterministic  large-scale  roughness  [8].  Thus,  zq  (r)  effectively  detrends  the  data 
and  ensures  (z)  =  0.  The  root-mean-square  (RMS)  surface  roughness  (sometimes  termed  interface  width) 
( 7  is  given  by  the  standard  deviation  of  fz(() 

(2) 

In  many  cases  the  height  PDF  of  real  random  rough  surfaces  is  Gaussian  in  form: 

/*(o=Rexp(“s?)  ■  w 

This  is  a  consequence  of  the  central- limit  theorem,  which  predicts  that  the  PDF  /s(C)  °f  the  sum  of  n 

n 

independent  identically  distributed  random  variables,  S  m  X)  A^,  tends  toward  Gaussian  form  as  n  — ►  oo, 

i—1 

provided  that  A i  have  finite  mean  and  finite  variance. 

Surfaces  having  identical  PDFs  are  distinguished  from  one  another  by  higher-order  statistical  descriptors, 
which  may  be  expressed  in  either  spatial  or  spectral  domains.  In  particular,  Gaussian  surfaces  (those  drawn 
from  an  ensemble  having  entirely  Gaussian  joint  PDFs)  are  fully  described  by  second-order  statistics.  In 
the  spatial  domain  second-order  statistics  are  expressed  by  the  auto  covariance  function  Cq: 

Co  (ri,  r2)  =  (z  (n)  z  (r2))  €  [-a2,  a2]  .  (4) 

Equivalently,  spatial  second-order  statistics  can  be  expressed  by  the  autocorrelation  function  C .  which  is 


3 


simply  the  autocovariance  function  normalized  by  the  variance  of  fz(() 


C(  ri,r2) 


{z  (rQ  z  (r2)) 

n-2 


€[-1.1]. 


(5) 


For  homogeneous,  isotropic  random  processes  C(i*i,r2)  is  a  function  only  of  the  distance  between  the 
two  points  C  ( R ),  where  R  =  |ri  —  T2  \  is  termed  the  translation  (sometimes  termed  the  lag  or  slip).  For  such 
surfaces,  a  correlation  length  £  can  be  defined  as  the  translation  at  which  C  (. R )  first  falls  from  its  maximum 
C  (0)  =  1  to  a  specified  value.  Here,  £  is  defined  as  the  translation  for  which  C  {£)  —  1/e.  In  other  words, 
the  correlation  length  is  a  measure  of  the  width  of  the  central  peak  of  the  autocorrelation  function.  While 
a  truly  random  process  will  have  some  value  of  translation  for  which  C  (R)  decreases  to  zero,  C  (R)  need 
not  be  a  strictly  decreasing  function  of  translation.  In  fact,  the  surfaces  of  interest  here  are  particularly 
distinguished  by  the  existence  of  nonzero  correlations  for  large  values  of  translation  (as  discussed  in  Sec. 
II C). 

The  structure  function  (also  termed  the  Allan  deviation,  the  height-height  correlation  function,  the 
mean-square  increment  function,  or  the  variogram)  is  a  spatial-domain  second-order  statistic  equivalent  to 
the  autocorrelation  function 

S  (R)  =  ([z  (r)  -z(  r  +  #)]2)  ,  (6) 

and  is  related  to  the  autocorrelation  function  by 


S  (R)  =  2a2  [1  -C  (i?)]  . 


(7) 


In  the  spectral  domain,  second-order  statistics  are  given  by  the  PSD.  For  homogeneous,  one-dimensional 
random  processes,  the  linear  PSD  p  (kx)  is  related  to  the  autocorrelation  function  and  vice  versa  by  a  Fourier 
transform  through  the  Wiener-Khint chine  relation: 

2  poo 

P  ( kx )  =  —  J  C(R)e~lkxRdR  ,  (8a) 


/OO 

P  (K), 

-oo 


llkxRdkx 


(8b) 


where  kx  £  (—00,  00)  is  the  one-dimensional  spatial  wavenumber.  When  expressed  in  radial  coordinates  for 
isotropic  two-dimensional  random  processes,  the  radial  PSD  P  (k)  is  related  to  the  autocorrelation  function 
and  vice  versa  by  a  Hankel  transform 


1  r°° 

P(k)  =  —  C  (R)  J0  (kR)  RdR  , 
J 0 

/»oo 

C{R)  =  2tt  /  P  (k)  J0  ( kR )  kdk  , 
Jo 


(9a) 


(9b) 


where  Jo  is  the  zeroth-order  Bessel  function  and  k  £  (0,  00)  is  the  radial  spatial  wavenumber  k  =  ||k||  = 

y/kl  +  fcg. 


B.  Power-Law  Rough  Surfaces 

Spectral  estimates  made  from  measurements  of  sea-floor  topography  are  often  well  modeled  by  PSDs  that 
are  inverse  power-law  functions  of  radial  spatial  wavenumber  P  (k)  oc  k~ 72  with  72  >  2  [9,  10,  11,  12,  13,  14] 
and  phase  spectra  that  are  uniformly  distributed  on  (0,  2tt]  [13,  14].  Because  of  the  tractability  of  scattering 
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theories  based  on  such  spectra,  characterizing  sea-floor  topography  parametrically  based  on  an  inverse- 
power-law  PSD  is  the  dominant  paradigm  in  underwater  geoacoustics,  among  other  areas  of  study. 

A  pure  power-law  PSD,  for  an  isotropic  surface,  is  given  by 


P  ^  =  (h0k)y2  f°r  k  G  ^°’  °°^  ’  (10) 

where  ho  is  a  reference  length,  usually  taken  to  be  unity,  w 2  is  the  spectral  strength  in  units  of  (length)4 
(i.e.,  the  value  of  the  PSD  at  k  =  h q1),  and  72  £  (2,  00)  is  the  two-dimensional  spectral  exponent. 

In  one  dimension,  (i.e.,  for  a  profile  measured  along  a  surface)  the  linear  PSD  is  given  by 

p(kx)  =  (hoZir  for  fc*€(~00’00)  ’  (11) 


where  71  =  72  —  1,  and 


w  1  = 


Vll  (hO 

h0  r  (i72) 


W  2 


as  shown  in  Refs.  15  (pp.  12-15)  and  2  (the  Appendix  A,  pp.  2898-2899).  (Note  that  both  of  the  cited 
references  contain  typographical  errors;  the  results  reported  here  are  correct.)  Regardless  of  functional  form, 
the  radial  PSD,  P(k),  and  linear  PSD,  p  (kx),  of  an  isotropic  stochastic  rough  surface  are  related  by  the 
Abel  transform  [see  Refs.  10  (the  Appendix)  and  15  (pp.  7-12)]. 

Note  that  p  (kx)  is  truly  a  power-spectral  density  in  that  it  has  units  (length)3,  which  are  units  of  “power” , 
(length)2,  per  elemental  unit  of  wavenumber,  (length)-1.  Likewise,  P  (k)  has  units  (length)4,  which  are  units 
of  “power”,  (length)2,  per  elemental  unit  wavenumber,  (length)-2.  Consequently,  for  profiles, 


a 


2 


(12) 


and,  for  surfaces, 

00 

a2  =  j  p  (kx)dkx  ,  (13) 

—  OO 

where  cr2  corresponds  to  the  total  “power”  of  z(r). 

General  bounds  on  spectral  exponent  are  expressed  in  terms  of  c),  the  topological  dimension,  by  7^  £ 
(D,oo),  where  D  =  1  or  2.  Topological  dimension  is  related  to  d,  the  dimensionality  of  the  embedding 
Euclidean  space,  by  d  =  D  +  l,  such  that  the  mapping  describing  the  topography  {z  (r)  :  r  £  IP}  corresponds 
to  the  graph  of  the  mapping  {(r,  z  (r))  :  r  £  M5},  which  is  in  M.d. 

It  is  nonphysical  for  the  PSD  of  surface  topography  to  follow  a  power-law  form  for  all  k.  The  resulting 
surface  has  infinite  RMS  roughness  because  of  the  singularity  of  the  PSD  as  k  — >  0.  Likewise,  this  nonphysical 
nature  is  reflected  in  the  autocorrelation  function  [see  Eq.  (15)  of  Ref.  2],  which  results  from  the  divergence 
of  the  integrals  relating  C(R)  to  the  PSD  of  power-law  form.  (However,  the  structure  function  [Eq.(6)] 
is  defined  for  surfaces  having  a  power-law  PSD  [16].)  Physically  realizable  random  rough  surfaces  must 
therefore  have  PSDs  that  depart  from  a  simple  power-law  form  for  small  k.  A  variety  of  PSDs  for  which 
power-law  behavior  has  a  low-frequency  cutoff  or  is  limited  to  finite  bandwidth  have  been  proposed;  some 
of  which  are  described  below  and  plotted  in  Fig.  1.  Particular  attention  has  been  given  to  PSDs  for  which 
C(R),  as  given  by  Eqs.(8  and  9),  can  be  given  an  analytical  expression.  This  is  of  interest  because  it  is  C(R) 
rather  than  P(k)  that  enters  directly  into  small-slope  calculations  of  scattering  from  rough  surfaces  [2]. 

The  simplest  way  to  achieve  a  PSD  having  the  needed  small- fc  cutoff  of  power-law  behavior  is  restricting 
the  nonzero  portions  of  a  power-law  PSD  to  bandwidth  [&_,  fc+]  where  k-  >  0  and  k+  may  be  infinite.  The 
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FIG.  1:  Power-spectral  density  as  a  function  of  radial  spatial  wavenumber  for  each  of  the  forms  described  in  Sec.  II B  are 
compared  for  72  =  3.0  and  u>2  =  2  •  10-5  m4.  Small-fc  cutoff  parameters  have  been  adjusted  to  yield  approximately  equal 

cutoff  points  for  each  PSD. 


resulting  spectrum  with  ‘box-car’  cutoff  is  given  by 


p{k)  = 


W  2 

(h0kyi 

0 


k-  ^  k  ^  fc+ 
otherwise 


(14) 


and  is  depicted  in  Fig.  1.  By  forcing  PSD  to  zero  at  k  —  0,  Eq.(14)  assures  that  the  resulting  surface  has 
no  constant  offset  in  height. 

As  an  alternative  to  the  abrupt  transitions  of  the  PSD  proposed  in  Eq.(14),  a  form  having  finite  power  at 
small  k,  which  may  better  represent  PSD  observed  both  in  shallow-water  (Ref.  17,  p.  365)  and  deep-water 
ocean  bottoms  [10,  18],  can  be  obtained  from  a  power-law  spectrum  with  an  ‘algebraic’  cutoff  at  small  k 
[2,  10]: 


P(k)  = 


w2 

r  9  9-172/2  5 

[( h0k )2  +  (h0kL)2\ 


(15) 


as  shown  in  Fig.  1.  A  power-spectral  density  of  a  similar  form  to  Eq.  15  can  be  approximated  by  a  piece- wise 
continuous  function  (not  shown  in  Fig.  1)  and  generalized  to  include  a  large- wavenumber  cutoff  [19],  which 
is  helpful  for  modeling  real  PSDs  measured  on  a  regular  array  of  points  (see  Sec.  Ill  A  1): 


p(k)  = 


W  2 

hok'l2 
w  2 

hok'Yz 

0 


k  <  kL 
kL  ^  k  ^  kH 
k  >  ku 


(16) 


Surfaces  described  by  Eqs.(15  and  16)  have  nonzero  mean  height,  which  can  be  removed  by  detrending  as 
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described  in  Sec.  IIIC  1. 

Observed  PSDs  sometimes  display  positive  slope  (power  decreasing  with  wavenumber)  below  a  cutoff 
frequency,  which  may  be  an  artifact  of  detrending  the  measured  profiles  [20] .  Power-spectral  density  of  this 
form,  shown  in  Fig.  1,  can  be  realized  by  a  power-law  spectrum  with  a  ‘Gaussian’  cutoff  at  small  k ,  which 
results  from  subtracting  the  Gaussian-weighted  moving  average  of  the  topography  from  a  surface  having 
power-law  PSD,  as  described  in  Ref.  21: 

p(fc)  =  (^{1"exp["(fcai)2/2]}2  •  (17) 

However,  C(R)  corresponding  to  this  PSD  is  not  analytically  expressible. 

A  similar  low-frequency  cutoff,  but  with  C(R)  that  is  analytically  expressible,  is  achieved  by  a  PSD 
written  as  the  difference  between  two  power-law  PSDs  having  algebraic  small- A:  cutoffs  [22,  23]. 

P(k)  =  P1(k)-P2(k)  ,  (18a) 


Pi  (k)  = 


h 


1  +  (fy^O^) 


72/2 


P2  (k)  = 


bi 


1  + 


[  b2hpk\ 

V  ) 


2"|  72/2  ’ 


(18b) 

(18c) 


as  shown  in  Fig.  1.  The  spectral  strength  w 2  and  the  location  of  the  spectral  peak  kp  are  related  to  the 
parameters  of  the  difference  spectrum  0,2 :  f>i,  and  62  through 


61 


w2&22 
1  -  a?  ’ 


(19a) 


U  1  ay^-1 

^0  kp  l  ^  _  ^1+72/2  ^ 

where  a 2  determines  the  steepness  of  P(k)  about  k  =  kp. 

While  spectral  exponent  72  is  often  assumed  constant  throughout  the  region  of  support  of  P(k),  it  has 
been  found  that  72  of  real  surfaces  sometimes  varies  with  k  [24,  25,  26].  Herzfeld  et  al.  have  observed  that 
spectral  exponent  is  greater  at  large  k  than  at  small  k  for  measurements  of  deep-water  ocean  bottoms  made 
over  length  scales  5  m-80  km  and  1  m-600  km  (see  Ref.  25  and  Refs,  therein).  At  much  smaller  scales, 
0.05  mm-3.5  m,  Briggs  et  al.  have  found  similar  results  for  measurements  of  shallow-water  granular  ocean 
bottoms.  This  increase  in  spectral  slope  at  large  k  is  often  associated  with  grain  size  of  the  constituent 
materials  [24,  27].  Measured  values  of  w<i  and  72  for  real  ocean  bottoms  are  described  in  Sec.  HE. 


C.  Fractal  Rough  Surfaces 

That  real  topography  is  taken  to  be  the  realization  of  a  random  process  having  P  (k)  oc  k~ 72  has 
significant  implications  for  how  it  is  characterized.  Power  spectra  of  this  form  are  frequently  encountered  in 
time-series  analysis,  where  random  sequences  z(x),  (D  =  1),  having  spectra  p  (kx)  oc  k~lx ,  71  £  (0,  00)  are 
termed  1  /fa  or  ‘fractional’  noise. 

The  range  of  random  sequences  of  this  form  is  bound  by  Gaussian  white  noise,  having  ji  =  0,  and 
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Brownian  motion  (realizations  of  the  Wiener  process,  sometimes  termed  Brown  noise),  having  71  =  2.  These 
two  bounding  sequences  are  related  in  that  Gaussian  white  noise  is  (formally,  over  smoothed  intervals)  the 
derivative  (or  ‘increments  process’)  of  Brownian  motion. 

Brownian  motion  (and  consequently  Gaussian  white  noise)  can  be  defined  axiomatically  [6,  28]. 

Theorem  1.  The  process  {z}  is  a  continuous  Brownian  motion  if  V  Ax  the  increments  A z(x)  = 
z  (pc  +  Ax)  —  z  (x)  satisfy 

1.  Gaussian  distribution 

2.  zero  mean 

3.  variance  oc  Ax 

4-  successive  increments  A z  (x)  and  A z  (x  +  Ax)  are  uncorrelated  .[84] 

An  axiomatic  definition  for  all  1  /  fa  noise  is  arrived  at  by  generalizing  the  axiomatic  definition  of  Brow¬ 
nian  motion  [29]. 

Theorem  2.  The  process  {z}  is  a  continuous  fractional  Brownian  motion  (FBM)  if\/  Ax  the  increments 
A z  (x)  =  z  (x  +  Ax)  —  z  (x)  satisfy 

1.  Gaussian  distribution 

2.  zero  mean 

3.  variance  oc  Ax2^;  where  S)  G  (0, 1]  is  the  Hurst  exponent 

4-  successive  increments  A z  (x)  and  Az  (x  +  Ax)  are  correlated  with  correlation  coefficient  given  by 
22^  =  2  +  2£  where  f  G  (-§,  l)  . 

Such  random  processes  are  termed  fractional  Brownian  motion  (FBM)  because  of  their  noninteger  spec¬ 
tral  exponents.  The  spectral  exponent  71  of  FBM  can  be  expressed  in  terms  of  the  Hurst  exponent  by 
71  =  2$)  +  1  (Ref.  30,  p.  67).  The  derivative  (increments  process)  of  FBM  is  termed  fractional  Gaussian 
noise  (FGN).  Together  FBM  and  FGN  span  an  important  segement  of  all  1  / fa  noise,  with  FBM  having 
71  G  (1,3]  and  FGN  having  71  G  [—1, 1].  The  only  Gaussian  process  that  is  both  self-affine  and  has  sta¬ 
tionary  increments  is  FBM  [31].  Fractional  Brownian  motion  (though  not  FGN)  can  be  generalized  from 
random  sequences  (D  =  1)  to  random  fields  (D  =  2)  and,  thus,  can  serve  as  a  model  for  the  topography  of 
random  rough  surfaces  having  a  power-law  PSD. 

It  is  particularly  notable  that  FBM  lacks  an  inherent  length  scale.  Stochastic  processes  that  have  this 
property  are  termed  self-affine  and  are  characterized  by  their  Hurst  exponent  S)  G  (0, 1].  A  real- valued 
stochastic  process  is  self  affine  with  Hurst  exponent  S)  if  its  increments  process  Az  (d)  =  z  (r  +  d)  —  z  (r) 
satisfies  the  homogeneity  relation 


{Az  (Xd)}  =  {X*Az  (tf)}  V  A,  r,  and  d  ,  (20) 

where  A  is  a  constant  scale  factor,  and  =  denotes  statistical  equality,  or,  more  exactly,  equality  of  the 
finite-dimensional  distributions  [32,  33].  This  homogeneity  relation  implies  that  empirical  histograms  of 
increments  processes  should  be  nearly  equal  when  scaled  according  to  Eq.  20.  Because  Eq.  20  is  valid  for 
all  d,  the  homogeneity  relation  also  holds  for  the  process  itself  [31] 


(z(Ar)}  =  {A*z(r)}  VA. 


(21) 


Consequently,  self-affine  surfaces  are  (statistically)  invariant  under  transformations  of  the  form 


(x,y,z)  =>  (A*,  Xy,  X^z)  . 


(22) 


If  =  1,  FBM  is  self-similar  and  appears  (statistically)  identical  when  observed  at  all  length  scales.  This 
is  impossible  for  single- valued  topography  that  can  be  expressed  by  the  mapping  described  in  Sec.  II  A,  for 
which  smaller  features  implictly  have  steeper  slope  than  larger  features  [34,  35].  If  S )  <  i,  FBM  appears 
(statistically)  identical  (e.g.  self-similar)  only  if  the  z- axis  is  scaled  independently  of  the  other  axes. 

A  property  of  self- affinity  is  dependence  of  surface  measures  (or,  more  exactly,  the  expectation  values  of 
these  measures  for  the  process)  on  observation  length,  i.e.  FBM  is  nonstationary.  For  example,  the  RMS 
roughness  [32,  33]  and  the  structure  function  obey  scaling  relationships 


a(L)  =  a(L0) 


ft 


(23a) 


S(L)  =  S(L0)^y  ,  (23b) 

where  L  is  the  length  of  an  interval  [0,  L\  contained  in  the  total  interval  [0,  Lq\  of  length  Lq.  Only  in  the  case 
that  S)  — >  0  are  surface  measures  independent  of  profile  length  and  is  FBM  a  stationary  random  process. 

In  the  case  of  FBM,  £ j  can  be  thought  of  as  a  parameter  that  determines  the  “rate  at  which  roughness 
changes  scale,”  [36]  such  that  smaller  values  of  S)  correspond  to  larger  variations  in  small-scale  topography. 
This  behavior  can  be  categorized  via  the  concept  of  ‘persistence’,  which  describes  the  extent  to  which 
consecutive  elements  of  increments  process  of  z  (r)  are  correlated.  Realizations  of  FBM  with  D  =  1  have 
S)  G  (0, 1).  For  S)  G  (0, 1/2)  FBM  is  antipersistent,  with  negative  correlation  between  consecutive  increments 
and  noise-like  spatial  variation.  For  S)  G  (1/2,1)  FBM  is  persistent  with  positive  correlation  between 
consecutive  increments  and  more  regular  spatial  variation.  In  the  case  that  Sj  ='  1/2,  FBM  reduces  to 
conventional  (or  ‘chaotic’)  Brownian  motion  with  uncorrelated  increments. 

The  graphs  of  realizations  of  FBM  are  one  of  a  class  of  geometries  termed  fractals — a  term  first  coined 
by  Mandelbrot  [28]  to  describe  profiles  and  surfaces  that  display  complexity  not  well  described  by  their 
topological  or  Euclidean  dimensionality.  Fractals  can  be  formally  (but  not  inclusively)  defined  as  sets  for 
which  the  Hausdorff-Besicovitch  dimension  H  exceeds  the  topological  dimension  T)  (Ref.  6,  p.  15).  Fractals 
may  be  either  deterministic  or  stochastic.  When  encountered  in  nature,  fractals  are  almost  exclusively 
stochastic.  In  other  words,  they  are  realizations  of  random  processes,  like  FBM. 

The  Hausdorff-Besicovitch  dimension,  H  G  [h,d),  is  a  measure  of  the  space- filling  quality  of  a  fractal 
object.  A  smooth  profile  (D  =  1,  d— 2)  has  H  —  1  and  a  smooth  plane  (D  =  2,  d— 3)  has  H  —  2.  While 
computing  H  of  a  set  is  complicated  (Ref.  6,  pp. 362-365, 458-459),  H  of  self-affine  profiles  and  surfaces 
is  equivalent  to  more  simply  computed  dimensions,  such  as  the  box  dimension  [37].  For  0—1,  the  box 
dimension  B  is  given  by 


—  logn 


B  —  lim 

x^o+  logx 


(24) 


where  n  is  the  number  of  boxes  of  linear  size  x  needed  to  cover  the  set  of  points  S  G  M  defining  the  profile. 

The  Hausdorff-Besicovitch  dimension  is  essentially  a  local  property  of  a  profile  or  surface  [Refs.  6  (p. 
373)  and  38].  In  contrast,  the  Hurst  exponent  characterizes  long-range  correlation  on  a  global  scale.  While 
not  true  for  all  fractals  [38],  for  self-affine  fractals,  such  as  FBM,  the  two  aspects  are  interrelated.  For 
self- affine  random  fractals,  the  Hurst  exponent  (global  property)  is  related  to  the  Hausdorff-Besicovitch 
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dimension  (local  property)  by  [see  Refs.  39  (pp.  301-325)  and  6  (p.  373)] 

H  =  d-S)  (25) 

such  that  smaller  S)  corresponds  to  larger  H  and  a  rougher  surface. 

Topothesy  (or  ‘crossover  length’)  r  is  a  measure  of  the  ‘strength’  of  a  fractal  and  is  defined  as  the 
horizontal  distance  (i.e.,  R)  between  two  points  connected  by  a  line  of  RMS  slope  one  radian  [11,  16,  40]. 
For  an  isotropic  rough  surface,  r  is  identical  for  the  surface  (D  =  2)  and  profiles  formed  by  cutting  the 
surface  (D  =  1).  Topothesy  can  be  alternatively  defined  in  terms  of  the  structure  function 

S  (r)  =  t2  V  D  .  (26) 

The  RMS  roughness  measured  over  a  profile  of  length  L,  a(L),  is  related  to  topothesy  through  the  scaling 
argument 

*(L)=t1-*L*  ,  (27) 

where  a  (r)  =  r. 

Together  the  Hausdorff-Besicovitch  dimension  (or,  alternatively,  the  Hurst  exponent)  and  topothesy 
provide  a  complete  description  of  an  isotropic,  self- affine  stochastic  fractal  surface. 

1.  Spectra  of  fractals 

As  conjectured  by  Mandelbrot  and  Van  Ness  (Ref.  29,  p.  436)  and  later  shown  by  Flandrin  [41,  42],  the 

formally  averaged  PSD  of  FBM  is  given  by  the  power-law  P  (k)  oc 

The  Hausdorff-Besicovitch  dimension  of  FBM  (of  arbitrary  D)  is  given  by  [16,  43] 

H  —  1  +  U  +  (c)  —  y^ )  /2  for  y?>  E  (c) ,  U  +  2]  (28) 

Because  Eq.  (28)  breaks  down  for  y^,  >  D  +  2,  it  is  often  argued  that  y2  E  (D,D  +  2).  However,  y^,  >  D  +  2 

simply  corresponds  to  nonfractal  differentiable  profiles  or  surfaces  of  H  —  U  having  smoothness  increasing 
with  y 3  [24,  44].  Such  a  conclusion  is  further  confirmed  by  experimental  measurements  of  spectral  exponents 
that  are  larger  than  the  proposed  limits  (see,  e.g.,  Ref.  12). 

Profiles  and  surfaces  can  be  categorized  based  on  fractal  dimension,  as  determined  by  spectral  exponent. 
For  profiles  (D  =  1),  yi  E  (1,  3]  =>  H  E  (2, 1]: 

•  yi  =  1,  H  =  2:  extreme  fractal  (borderline  area  filling) 

•  yi  =  2,  H  =  1.5:  Brownian  fractal  (graph  of  one-dimensional  Wiener  process) 

•  yi  =  3,  H  =  1:  marginal  fractal  (borderline  nonfractal)  . 

Five  profiles  that  are  realizations  of  FBM  are  shown  in  Fig.  2.  Fractal  dimensions  of  these  profiles  range 
from  2  to  1  in  uniform  increments  of  0.25. 

For  surfaces  (D  =  2),  y2  E  (2, 4]  =>  H  E  (3,  2]: 

•  y2  =  2,  H  =  3:  extreme  fractal  (borderline  volume  filling) 

•  y2  =  3,  H  =  2.5:  Brownian  fractal  (graph  of  two-dimensional  Wiener  process) 

•  y2  =  4,  H  =  2:  marginal  fractal  (borderline  nonfractal)  . 
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FIG.  2:  Realizations  of  FBM  for  the  range  of  Hausdorff-Besicovitch  dimension  H  E  [2, 1].  Hausdorff-Besicovitch  dimensions 
are  2  (a),  1.75  (b),  1.5  (c),  1.25  (d),  1  (e).  The  y- axis  is  in  normalized  units, as  the  amplitudes  of  the  zero-mean  profiles  are 
normalized  by  cr.  The  a?- axis  is  in  arbitrary  units  (sample  number) 


The  topothesy  of  a  FBM  profile  (D  =  1)  or  surface  (D  =  2)  is  related  to  power-law  parameters  by  [45] 


W27rX)/22ar  (a/2) 
ho  (7z>  -^)r(7t)/2) 


for  7a  E  (D,D  +  2)  , 


where  a  =  D  +  2  —  73  such  that,  for  isotropic  surfaces  (D  =  2), 


w2 


7T(72-2)r(g)  2 

22-72 r  (^)  (2tt)72 


(29) 


(30) 


Band-limited  fractals 

The  power-law  PSD  of  self- affine  stochastic  fractal  profiles  and  surfaces,  such  as  FBM,  is  problematic 
when  the  power-law  functional  form  is  assumed  to  exist  for  all  k,  as  discussed  in  Sec.  II B.  FBM  cannot  be 
easily  characterized  using  conventional  metrics,  such  as  those  described  in  Sec.  II  A.  In  particular,  realiza¬ 
tions  of  FBM  possess  mulitscale  roughness  such  that  they  are  continuous  but  nondifferentiable.  Additionally, 
they  display  long-range  correlations  over  a  range  of  translations  limited  only  by  the  physical  extent  of  the 
particular  realization.  Self-affinity  violates  assumptions  of  (wide-sense)  stationarity  due  to  dependence  of 
measured  quantities  (e.g.,  cr)  on  sample  size,  while  long-range  correlations  violate  assumptions  of  ergodicity 
due  to  the  inability  of  any  finite  spatial-averaging  domain  to  be  large  enough  to  recover  all  information 
about  the  random  process  [46]. 

Real  surfaces  can  be  self-affine  only  over  a  “finite  hierarchy  of  scales”  [45].  Such  surfaces  can  be  said 
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to  be  “approximately  self-affine”  [45]  or  to  posses  “limited  self-affinity”  [47].  At  large  scales,  the  finite 
extent  provides  a  small- fc  limit.  At  small  scales,  the  molecular  or  granular  structure  provides  a  large- k 
limit.  Alternatively,  small-scale  limits  of  self- affinity  can  result  from  the  limited  resolution  of  manufacturing 
processes,  in  the  case  of  fabricated  rough  surfaces,  or  as  an  artifact  of  finite- resolution  measurement  methods, 
when  data  from  from  real  surfaces  are  analyzed. 

Rough  surfaces  are  therefore  better  modeled  as  approximately  self-affine,  having  band-limited  PSDs, 
such  as  those  given  in  Sec.  II B.  Doing  so  has  a  number  of  advantages.  In  particular,  conventional  metrics, 
such  as  RMS  roughness,  cr,  and  autocorrelation,  C  (R),  which  are  not  well  defined  for  FBM,  can  be  defined 
unambiguously  for  band-limited  FBM.  It  has  been  shown  that  both  FBM  profiles  [41]  and  surfaces  [48] 
having  power-law  PSD  restricted  to  a  finite  band  (with  some  additional  conditions  on  the  falloff  outside  of 
the  band)  are  wide-sense  stationary.  Likewise,  band-limited  FBM  profiles  and  surfaces  are  differentiable, 
due  to  the  large- A:  cutoff. 

However,  limiting  the  hierarchy  of  scales  of  self-affinity  for  a  self-affine  fractal  can  alter  the  Hausdorff- 
Besicovitch  dimension  H  from  that  of  the  equivalent,  ideally  self-affine  fractal  [49,  50,  51,  52,  53].  The 
relationship  between  H  and  7^  given  in  Eq.  (28)  holds  only  for  ideally  self- affine  fractals.  Band  limiting 
in  the  spatial- wavenumber  domain  alters  the  form  of  the  structure  function,  S  (R),  such  that  the  linear 
relationship  between  log  S'  (R)  and  logR  extends  only  over  the  scaling  interval  R  G  (R_,R+),  defined  by 
[50]: 


r  (h)  l1/(3  7l)  2 

_7i-ir(7l/2)J  G~ 


(31a) 


R+  = 


r(7i/2)  l1/(7l_1)  2 

_v/5FT  (tf-l)J  GG  ’ 


(31b) 


where  H  is  the  Hausdorff-Besicovitch  dimension  of  an  ideal  self-affine  fractal  having  the  same  spectral 
exponent  ( H  =  (5  —  7i)/2). 

Similarly,  topothesy,  as  defined  by  Eq.  (26),  is  a  function  of  the  shape  and  location  of  small-  and  larg e-k 
cutoffs,  which  alter  the  form  of  the  structure  function  [19].  The  small- fc  cutoff  is  manifest  as  oscillations  in 
S(R)  (or,  alternatively,  C(R))  at  large  R  [7,  50].  The  larg  e-k  cutoff  is  also  manifest  as  oscillations,  which 
are  noticeable  only  for  large  a  [50] . 

An  alternative  definition  of  r,  based  on  an  analytic  expression  for  S(R ),  is  independent  of  wavenumber 
cutoffs.  The  band  limited  power-law  PSD  given  by  Eq.  (14)  has  S(R)  that  is  expressible  as  an  asymptotic 
expansion  [50,  54]  the  first  term  of  which  is  S(R)  of  an  ideally  self- affine  fractal 


S{R)  =  rZ)+2-70  r^-v  . 


(32) 


Defining  r  by  the  above  expression  leads  to  a  form  that  is  independent  of  the  spectral  cutoffs. 


D.  Parameter  Study:  Realizations  of  Band-Limited  Fractal  Surfaces 

This  section  presents  a  series  of  realizations  of  a  band-limited  FBM  surface.  These  are  shown  in  Figs. 
3-6  and  described  in  Table  I.  These  figures  are  grouped  into  two  series,  as  indicated  in  Table  I:  Fig.  3 
is  a  movie  (with  frames  shown  individually  in  the  Appendix  as  Figs.  13-18)  illustrating  the  the  effect  of 
increasing  72;  Figs.  4-6  show  the  dependence  on  low- wavenumber  cutoff  [here  ka l  of  Eq.  (17)]. 

Varying  72  and  kaL  produces  changes  in  the  realized  surface  that  can  be  appreciated  visually  (the 
parameter  simply  acts  as  a  scale  factor).  In  the  movie  shown  in  Fig.  3  (and  also  Figs.  13-18  of  the 
Appendix),  72  is  varied  for  a  particular  realization  while  the  other  two  parameters  are  fixed.  Varying  72 
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TABLE  I:  Realizations  of  a  band- limited  fractal  surface.  All  cases  have  /iq  =  lm  and  W2  =  10  5  m4. 


72 

2.0 

2.4 

2.8 

3.2 

3.6 

4.0 

2.8 

2.8 

2.8 

kdL 

0.05 

0.05 

0.05 

0.05 

0.05 

0.05 

0.01 

0.05 

0.10 

Figure 

13° 

14a 

15a 

16° 

17a 

18° 

4 

5 

6 

°Collected  in  Fig.  3  as  one  frame  of  an  embedded  QuickTime  animation. 


introduces  a  change  in  the  overall  length  scale,  but  as  the  idealized  distribution  of  heights  is  Gaussian,  the 
range  of  heights  can  be  approximately  standardized  by  subtracting  the  mean  height  of  the  entire  realization 
and  dividing  by  the  standard  deviation  of  the  heights.  In  each  figure  (including  the  figures  in  the  appendix), 
the  x  and  y  axes  are  in  units  of  meters  and  the  z-axis  (represented  by  the  color  bar)  is  in  normalized  units 
arrived  at  by  dividing  by  a.  As  such,  realizations  are  presented  using  the  same  color  scale.  It  is  observed 
that  the  graininess  of  various  features  decreases  as  72  increases,  until  the  surface  looks  qualitatively  different 
from  the  72  =  2  surface.  This  is  related  to  the  behavior  of  the  correlation  function,  as  the  values  at  small 
lags  can  differ  considerably  for  smaller  values  of  72. 

The  particular  plotting  technique  chosen  results  in  certain  features  of  the  surface  being  perceived  as 
islands,  the  boundaries  of  which  have  a  level  of  detail  that  varies  with  the  fractal  dimension  of  the  surface. 
In  this  way,  the  presentation  is  similar  to  the  “zero  set”  of  the  elevation  data  formed  by  the  intersection  of 
the  surface  with  a  plane  oriented  parallel  to  the  nominal  surface  orientation.  The  fractal  dimension  of  the 
boundaries  of  the  islands  formed  is  thus  one  less  than  the  fractal  dimension  of  the  surface  itself.  [55] 

Figures  4-6  show  the  results  of  varying  the  low -k  cutoff.  In  comparison  with  Fig.  5,  which  has  the 
same  value  of  72  =  2.8  as  Figures  4  and  6,  there  is  an  observable  qualitative  effect  that  changing  the  cutoff 
produces  a  change  in  the  apparent  values  of  72,  the  high  cutoff  value  producing  a  realization  that  has  an 
apparent  72  greater  than  3. 

E.  Rough  Surface  Experiments  with  Physical  Models 

As  stated  in  the  Introduction,  rough-surface  experiments  with  physical  models  become  highly  useful 
when  a  systematic  field  experiment  is  impractical,  as  in  a  controlled  experiment  with  bottom  roughness  pa¬ 
rameters.  While  it  is  difficult  to  create  physical  models  that  are  exact  acoustic  analogs  of  real  ocean  bottoms 
over  a  complete  set  of  parameters  (such  as  density,  shear /compressional  speed,  and  shear/compressional  at¬ 
tenuation),  a  physical  model  can  be  constructed  to  exhibit  scattering  behavior  analogous  to  a  real  ocean 
bottom,  such  as  a  strong  dependence  on  a  particular  parameter. 

As  an  example,  Carolina  Coast  limestone  is  a  type  of  elastic  ocean  bottom  for  which  the  scattering 
behavior  has  a  large  predicted  dependence  on  72 .  At  two  sites  off  the  Carolina  Coast,  Soukup  and  Gragg 
[56]  extracted  values  of  w<i  and  72  using  a  geoacoustic  inversion  technique  that  analyzed  scattering-strength 
vs.  grazing- angle  curves  in  the  band  2-3.5  kHz  in  water  depths  of  230-310  m  under  the  assumption  that 
interface  scattering  wholly  dominated  any  volume-scattering  effects.  The  inverted  spectral  parameters  for  the 
Carolina  Coast  limestone  obtained  by  Soukup  and  Gragg  [56]  were  72  =  2.6  —  2.8  and  w<i  —  4- 10-4  —  9  •  10-4. 
The  scattering  strength  for  this  limestone  predicted  by  the  NRL  small-slope  model  [2]  is  plotted  in  Fig.  7 
as  a  function  of  grazing  angle  for  different  values  of  72  at  3500  Hz. 

Limestone  supports  both  compressional  and  shear  waves,  each  with  associated  critical  angles.  Particu¬ 
larly  characteristic  of  the  behavior  of  limestone  is  a  clear  stratification  of  scattering  strength  as  a  function 
of  72  below  the  lowest  (shear-speed)  critical  angle.  To  investigate  this  substantial  dependence  of  scatter¬ 
ing  strength  on  72,  the  predicted  behavior  of  PVC,  a  material  with  analogous  backscattering  behavior, 
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FIG.  3:  Embedded  QuickTime  animation  illustrating  the  changes  observed  in  realizations  of  band-limited  fractal  surfaces 
with  as  72  varies  from  2  to  4,  for  fixed  W2  =  1  •  10“ 5  m4  and  kai,  =  0.05.  (See  Table  I.) 


was  studied.  Below  the  compressional-speed  (p-wave)  critical  angle  for  PVC  at  300  kHz,  the  dependence 
of  scattering  strength  on  72  is  even  more  pronounced  than  for  limestone.  (Though  PVC  supports  shear 
waves,  the  speed  of  these  waves  is  such  that  no  shear- wave  critical  angle  exists.)  An  example  is  shown  in 
Fig.  8.  Scattering  strengths  from  a  rough  PVC  sheet  at  the  Allied  Geophysical  Laboratory /University  of 
Houston  were  obtained  and  the  mean  behavior  of  the  scattering  was  comparable  to  the  predictions  of  the 
NRL  small-slope  scattering  model  [57].  The  PVC  sheet  had  mean  values  of  W2  and  72  of  1.5  •  10-5  m4  and 
2.9,  respectively.  The  authors  are  currently  fabricating  rough  surfaces  from  PVC  that  will  represent  other 
values  of  and  72.  Future  research  will  also  introduce  different  materials  so  that  other  predicted  effects 
of  the  small-slope  theory,  such  as  local  minima  at  particular  grazing  angles,  can  be  verified.  While  the 
lowest-order  small-slope  predictions  shown  here  were  the  motivation  for  the  authors’  experimental  design, 
modeling  the  actual  scattering  behavior  for  such  surfaces  may  require  the  use  of  a  higher-order  model  for 
surface  scattering  and  a  model  of  scattering  in  the  subsurface  volume. 


III.  FABRICATION  OF  ROUGH  SURFACES 

Creating  a  physical  realization  of  a  random  rough  surface  is  a  two-step  process.  The  first  of  these  steps  is 
the  numerical  generation  of  the  topography.  Section  III  A  describes  how  this  can  be  achieved  by  a  spectral 
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FIG.  4:  Realization  of  a  band-limited  fractal  surface  with  W2  =  1  •  10  5  m4,  72  =  2.8,  and  kdL  =  0.01.  (See  Table  I.) 


method.  In  this  section  the  theoretical  effects  of  the  finite  size  of  the  realization  and  the  requirement  that 
the  surface  be  defined  only  over  a  discrete  grid  of  points  are  also  considered.  The  second  step  of  creating 
a  physical  realization  of  a  random  rough  surface  is  physically  manufacturing  the  generated  topography. 
Section  IIIB  describes  how  this  can  be  achieved  using  a  computer- numerically  controlled  milling  machine. 
In  practice,  the  topography  of  a  physically  manufactured  surface  must  be  measured  to  inspect  for  accuracy 
of  the  manufacturing  process.  Such  a  measurement  process,  using  a  touch-trigger  probe,  is  described  in  Sec. 

me. 


A.  Generating  Surface  Realizations 


When  modeled  as  a  random  field,  the  topography  of  a  particular  rough  surface  z  (r)  is  a  realization  of  a 
particular  two-dimensional  random  process.  Taking  FBM  as  a  model  of  random  rough  surface  topography, 
the  first  step  in  fabricating  a  rough  surface  is  generation  of  a  realization  of  FBM. 

The  earliest  methods  for  generating  realizations  of  FBM  utilized  the  Weierstrass-Mandelbrot  function 
[58].  This  extension  of  the  Weierstrass  function  [59,  60],  first  proposed  by  Mandelbrot  [28],  represents  profiles 
(0  =  1)  by 


K  [W  7)]  = 


^  [(1  -  eianx )  e^»] 

2^  a(2 -H)n 


(33) 


n=  —  00 
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FIG.  5:  Realization  of  a  band-limited  fractal  surface  with  W2  =  1  •  10_5  m4,  72  =  2.8,  and  kaL  =  0.05.  (See  Table  I.) 

where  a  >  1  is  a  parameter  and  (j)n  is  the  phase,  a  random  variable  uniformly  distributed  on  (0,  27t].  Equation 
(33)  describes  a  summation  of  sinusoids  having  the  geometric  spectrum  of  wavenumbers  kx  =  an,  the  so- 
called  Weierstrass  spectrum,  and  amplitudes  oc  kx  =  an.  The  Weierstrass-Mandelbrot  function  is  known  to 
be  a  fractal  [28,  58,  61]  and  has  a  mean  PSD  consistent  with  FBM  [58]. 

While  first  used  only  to  describe  fractals  generated  by  Eq.(33),  as  in  Ref.  58,  the  term  Weierstrass- 
Mandelbrot  fractal  has  been  informally  extended  to  describe  similar  summations  of  sinusoids  having  different 
wavenumber  spectra.  In  particular,  those  defined  on  spectra  of  algebraic  intervals,  such  as  the  Fourier 
spectrum  [47],  can  also  be  termed  Weierstrass-Mandelbrot  fractals.  Berry  and  Lewis  [58]  were  first  to  use 
the  Weierstrass-Mandelbrot  function,  defined  over  the  Weierstrass  spectrum,  as  models  of  FBM  (D  =  1). 
Equivalent  synthesis  over  the  Fourier  spectrum  is  discussed  in  Ref.  62  (pp.  76-77).  Kumar  and  Bodvarsson 
[63]  later  extended  these  methods  to  surfaces  (D  =  2),  both  for  Weierstrass  and  Fourier  spectra.  (See  also 
Ref.  62,  pp.  174-187.) 

When  defined  over  the  Fourier  spectrum,  the  Weierstrass-Mandelbrot  function  is  essentially  equivalent 
to  an  inverse  Fourier  transform  of  a  k- space  (spatial-wavenumber  domain)  description  of  FBM  consisting 
of  a  power-law  PSD  and  a  random  phase  spectrum.  Here  the  Fourier  method  is  chosen  for  the  numerical 
generation  of  surface  realizations  because  discrete  Fourier  transforms  on  finite,  uniform  grids  can  be  efficiently 
calculated  by  the  fast  Fourier  transform  (FFT). 
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x(m) 

FIG.  6:  Realization  of  a  band-limited  fractal  surface  with  W2  =  1  •  10-5  m4,  72  =  2.8,  and  kaL  =  0.1.  (See  Table  I.) 

1.  Spatial  sampling  and  surfaces  of  finite  extent 

For  the  numerical  generation  of  surface  realizations,  z  (r)  is  expressed  by  a  finite  two-dimensional  array 
representing  surface  heights  sampled  on  a  uniform  grid.  The  relation  between  the  spatial-sampling  interval 
Ax  and  the  maximum  resolvable  linear  spatial  wavenumber  kx  —  max (\kx\)  in  each  Euclidean  dimension 
(equivalent  to  the  Nyquist  frequency),  is  given  by  the  Whittaker-Shannon-Kotel’nikov  (WSK)  sampling 
theorem  [64] 

kx  =  7T /  Ax  .  (34) 

For  isotropic  surfaces  sampled  on  a  uniform  grid  ( Ax  =  Ay ),  the  maximum  radial  spatial  wavenumber  is 
k  —  V%kx. 

The  relation  between  the  linear  extent  of  the  surface  (region  of  support)  X  and  the  minimum  resolvable 
linear  spatial  wavenumber  kx  =  min(|A:x|)  in  each  Euclidean  dimension  (that  for  which  one  wavelength 
corresponds  exactly  to  the  linear  extent  of  the  surface)  is 

kx  =  2tt/X  .  (35) 

Clearly,  kx  <  kx  can  influence  the  piston  (zeroth-order  polynomial  term),  tilt  (first-order  polynomial  term), 
and  curvature  (second-  and  higher-order  polynomial  terms)  of  the  surface,  though  detrending  strongly 
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FIG.  7:  Lowest-order  small-slope  model  predictions  of  the  backscattering  strength  vs.  grazing  angle  at  three  values  of  72  for 
an  incident  wave  of  3  kHz  scattered  from  a  limestone  bottom  of  the  Carolina  Coast,  with  a  W2  value,  density, 
shear /compressional  speeds,  and  shear /compressional  attenuations  inverted  from  at-sea  experimental  data  [56].  Note  the 
stratification  of  scattering  levels  as  a  function  of  72  predicted  for  angles  below  the  shear-speed  (s-wave)  critical  angle. 


grazing  angle  (degrees) 

FIG.  8:  Lowest-order  small-slope  model  predictions  of  the  backscattering  strength  vs.  grazing  angle  at  three  values  of  72  for 
an  incident  wave  of  300  kHz  scattered  from  the  Houston  PVC  sheet,  with  a  W2  value  obtained  from  profilometry  and  the 
density,  shear /compressional  speeds,  and  shear /compressional  attenuation  measured  from  a  sample  of  the  material.  Note  the 
wide  (in  comparison  with  Fig.  7)  disparity  in  the  scattering  levels  for  angles  below  the  compressional-speed  (p-wave)  critical 
angle.  Verifying  this  disparity  is  one  goal  of  the  tank  experiments  currently  underway  at  NRL. 
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suppresses  kx  <  kx  [20].  For  surfaces  realized  on  a  square  surface  of  support  ( X  =  Y),  the  minimum  radial 
spatial  wavenumber  is  k  —  V2kx.  The  observed  PSD  of  a  physically  realized  power-law  surface  can  be 
modeled  by  Eq.  (14)  where  kL  =  k  and  kn  =  k. 


2.  Fractal  properties 


Discrete,  finite-support  realizations  of  FBM  are  approximately  self-affine,  band-limited  fractals.  The 
length-scale  range  x_  to  x+  over  which  a  such  a  profile  or  surface  is  self- affine  is  determined  by  (1)  the 
spatial  sampling  interval  Ax,  which  sets  the  large-/c  limit  x_  and  (2)  the  extent  of  the  region  of  support  X, 
which  sets  the  small- fc  limit  x+.  Specifically, 


Ax 

~Y  ’ 

(36a) 

X 

2 

(36b) 

3.  Computational  procedure 


In  this  work,  approximately  self-affine,  band-limited  realizations  of  FBM  are  generated  by  a  spectral 
method  in  which  a  representation  of  the  surface  Z( k)  is  created  in  k- space  (the  spatial- wavenumber  domain) 
and  then  transformed  into  a  representation  z( r)  in  real  space  (the  spatial  domain)  using  the  discrete  Fourier 
transform  (DFT). 

In  two  dimensions  the  DFT  is  defined  by  the  pair  of  equations 


Z  (kx ,  ky) 


Y  Y  z{nx,ny)e  ‘l(2n/N* 3kxnxe  i{2ir/Ny)kyny  ^ 


0^kx^Nx-l,0^ky^Ny-l 

otherwise 

(37a) 


Z  fax  i  ^y) 


Nx—1  Ny-1 

Y  E  Z  fax,  ky)  eK^/Nx)kxnXei{2,/Ny)kyny^ 
kx  =0  ky=0 


o, 


0  ^  nx  <  Nx  -  1, 0  <  ny  <  Ny  -  1 
otherwise 


(37b) 


for  a  finite-extent  first-quadrant  support  sequence  [65] . 

The  surface  z  (r)  must  be  purely  real,  which  requires  that  Z  (k)  be  conjugate  symmetric  [65] 

real  z  (nx,nv)  Z  (kx,ky)  =  Z*  (~kx,  -kx)  RnxxNv  (kx,ky)  , 

where  Z  indicates  here  that  the  function  is  periodic  with  a  period  of  Nx  X  Ny, 


X  ( kX  ,  ky) 


X  (kX  +  Nx,  ky)  —  X  ( kx ,  ky  +  Ny)  \/{Nx,  Ny)  , 


and  RnwxNv  (kx,ky)  is  defined  by 


Rnx  x  Ny  )kx ,  ky) 


1,  0  ^  kx  ^  Nx,  0  ^  ky  ^  Ny 
0,  otherwise 


(38) 


(39) 


(40) 
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In  general, 


Z(k)  =  |Z(k)|  e2niarsZW  =  sr  [Z{\z)\  +  9  [Z{ k)]  (41) 

where  P( k)  =  |Z(k)|2.  Here,  z( r)  is  isotropic  such  that  Zlk)  is  a  function  of  /.:  =  only.  The 

square  magnitude  of  Z( k)  is  chosen  to  follow  a  modified  power-law  form,  as  given  by  Sec.  II B.  The  phase  is 
chosen  to  be  random,  though  must  be  specified  such  that  it  meets  the  symmetry  requirements  of  Eq.  (38). 
Defining  a  as  an  antisymmetric  array  of  random  numbers  uniformly  distributed  on  (0, 1)  (as  described  in 
the  Appendix  B)  yields  a  simple  algorithm  for  generating  a  realization  of  the  surface: 

z( r)  =  T  \yp  (fc)e2™a]  ,  (42) 

where  T  denotes  the  Fourier  transform  operator.  While  this  algorithm  ensures,  by  construction,  that  the 
expected  value  of  the  PSD  of  the  generated  surface  is  P(k\  it  does  not  explicitly  ensure  that  the  resulting 
rough  surface  will  be  a  realization  of  a  Gaussian  random  process. 

A  similar  approach  to  generate  realizations  of  band-limited  FBM,  which  was  followed  by  Thorsos  [66], 
produces  realizations  for  which  the  expected  value  of  the  PSD  is  P(k),  and  which  also  have  surface-height 
distributions  that  are  guaranteed  to  be  Gaussian.  In  this  method,  z( r)  is  generated  by  the  algorithm 


*(r)  =  ? 


VpW) 


/jL  +  irj 


V2 


(43) 


where  fi  and  r]  are  uncorrelated  arrays  of  random  numbers  normally  distributed  with  mean  of  zero  and 
standard  deviation  of  one. 

An  equivalent  algorithm  to  that  given  in  Eq.  (43),  which  was  used  here,  can  be  obtained  by  modifying 
Eq.  (42)  such  that  the  magnitude  is  randomly  distributed  about  P(k)  according  to  a  Rayleigh  distribution 

z( r)  =  F  [y-  In  (/?)  P  (fc)e2,riaJ  ,  (44) 

where  f3  is  a  symmetric  array  of  random  numbers  uniformly  distributed  on  (0, 1)  (as  described  in  the 
Appendix  B). 

The  Fourier  transform  used  by  these  algorithms  is  a  one-to-one  linear  transformation  (though  it  is 
not  shift-invariant).  It  can  be  shown  that  a  linear  transformation  of  a  realization  of  a  Gaussian  random 
process  yields  a  sequence  that  is  also  a  realization  of  a  Gaussian  random  process  (Ref.  67,  pp.  340- 
341).  Consequently,  a  random  process  that  is  Gaussian  distributed  must  have  a  Fourier  transform  (k- space 
representation)  that  is  also  Gaussian  distributed  in  both  its  real  and  imaginary  components  [68].  The 
algorithm  given  in  Eqs.  (43)  and  (44)  therefore  yields  rough  surfaces  that,  by  construction,  are  realizations 
of  a  Gaussian  random  process.  In  contrast,  the  algorithm  given  in  Eq.  (42),  yields  rough  surfaces  that 
are  approximately  realizations  of  a  Gaussian  random  process  through  the  Central  Limit  Theorem.  The 
inverse  DFT,  Eq.  (37b),  of  a  random  complex- valued  random  process  in  k  space  represents  each  sample  of 
the  transformed  process  in  real  space  as  a  linear  combination  of  random  variables,  that  is,  the  individual 
samples  of  the  &-space  representation.  Provided  that  these  random  variables  are  mutually  independent  with 
zero  means  and  finite  variances  (Ref.  67,  pp.  225-230),  the  Central  Limit  Theorem  predicts  that  each 
sample  in  real  space  will  tend  toward  a  Gaussian  distribution  as  the  number  of  frequency  terms  in  the 
inverse  DFT  becomes  large.  While  the  Central  Limit  Theorem  is  often  approximately  valid  for  surprisingly 
small  numbers  of  random  variables,  the  breakdown  of  the  model  can  be  aggravated  by  the  power-law  form 
of  the  PSD,  which,  by  decreasing  amplitude  at  large  k  effectively  eliminates  the  influence  of  large- k  terms. 
By  using  the  algorithm  given  in  Eqs.  (43)  and  (44),  the  realization  effectively  simulates  a  larger  surface. 
However,  it  can  be  argued  that  Eq.  (42)  is  more  representative  of  the  central-limit  processes  through  which 
physical  surfaces  arrive  at  a  Gaussian  height  distribution. 


20 


B.  Manufacturing  Rough  Surfaces 


Various  approaches  have  been  followed  for  the  manufacture  scale-model  rough  surfaces  for  use  in  acousti¬ 
cal  experiments.  There  is  likewise  a  diversity  of  sources  of  fabricated  topography,  some  electing  to  reproduce, 
at  scale,  topography  measured  from  actual  surfaces  and  others  choosing  to  generate  arbitrary  surfaces  sat¬ 
isfying  certain  statistical  conditions. 

In  the  first  reported  work  of  this  type,  Proud  et  al.  [69]  built  an  approximately  Gaussian  rough  surface 
from  layers  of  cork  strip  fixed  to  plywood  substrate.  Later  Horton  et  al.  [70]  manufactured  surfaces  from 
low-density  expanded  polystyrene  based  on  aeromagnetic  contours  of  the  Canadian  Shield.  Though  not 
detailed,  such  a  machining  process  must  have  been  quite  difficult  given  the  state  of  the  art  at  the  time 
(1967).  Much  more  recently,  Gautier  and  Gibert  [36]  used  resin  printing  to  create  a  replica  of  a  natural 
fracture  surface  of  granite,  the  topography  of  which  was  measured  using  laser  profilometry. 

A  particularly  promising  technique,  which  was  adopted  for  the  work  described  here,  is  use  of  computer- 
aided  manufacturing  (CAM)  to  shape  surfaces  of  near  arbitrary  topography  from  blocks  of  raw  material.  Use 
of  this  fabrication  technique  for  acoustical  scale  models  was  first  described  by  Mellema  (Ref.  71,  pp.  33-35), 
but  was  used  earlier  to  fabricate  rough  surfaces  for  other  purposes  [72,  73].  Computer-aided  manufacturing 
uses  software  to  develop  fabrication  strategies  directly  from  computer  representations  of  a  desired  finished 
object.  Such  machining  strategies  take  the  form  of  a  series  of  tool  paths  that  describe  how  a  particular 
sequence  of  cutting  tools  must  be  moved  in  order  to  achieve  a  prescribed  result.  Tool  paths  are  expressed  in 
a  scripting  language  (e.g.,  G  code)  and  provide  the  instruction  sets  used  by  computer- numerically  controlled 
(CNC)  tools,  such  as  milling  machines. 

Milling  machines  are  described  by  the  number  of  independent  axes  through  which  they  can  move  a  cutting 
tool.  The  number  of  axes  determines  the  maximum  Euclidean  dimension  (d)  of  the  space  embedding  an 
object  that  such  a  mill  is  able  to  fabricate.  Manufacture  of  random-rough  surface  topography  (D  =  2,  d  =  3) 
as  described  in  Sec.  II A  requires  a  3-axis  mill  able  to  independently  move  along  each  of  the  Cartesian 
axes.  Additional  axes,  such  as  employed  in  5-axis  mills  [74],  are  required  only  to  fabricate  objects  of  higher 
topological  dimension  or  surfaces  that  are  not  single- valued  mappings. 

Machining  of  profiled  surfaces  is  typically  performed  with  a  ball- nose  end  mill,  as  opposed  to  the  more 
common  flat  end  mill.  Complete  fidelity  in  milling  topography  z  (r)  requires  that  w,  the  radius  of  the 
ball- nose  end  mill,  be  smaller  than  p,  the  minimum  positive  radius  of  curvature  of  the  surface.  At  any 
given  point  r,  the  radius  of  curvature  p  (r)  is  found  by  first  calculating  the  mean  curvature  kh  (r)  and  the 
Gaussian  curvature  kk  (r)  via 
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where  subscripts  in  these  expressions  indicate  partial  derivatives  with  respect  to  the  coordinate  axes  x  and 
y ,  evaluated  at  r.  In  terms  of  these,  (r)  =  kh  (r)  ±  y/ (r)  —  Kk  (r)  are  the  two  principal  curvatures 
of  z  at  r.  The  minimum  positive  radius  of  curvature  at  r  is  p  (r)  =  1  /«+  (r),  such  that  p=\  j  max  [ft+  (r)]. 

While  easily  expressed,  the  above  requirements  pose  the  problem  of  differentiating  a  fractal  surface  that, 
by  definition,  is  nondifferentiable.  This  problem  is  removed  in  practice  because  the  fractal  surface  is  given 
a  discrete,  finite-support  (and  therefore  band-limited)  representation  that  is  only  approximately  self-affine. 
Consequently,  derivatives  are  defined  [19,  24]. 

The  requirement  of  fidelity  limits  the  maximum  radius  of  the  tool  used  to  mill  the  surface,  but  the  nec- 
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FIG.  9:  Simulation  of  a  ball-nose  end  mill  tracing  the  evolute  tool  path. 


essarily  small  cutting-edge  length  of  such  a  tool  may  be  inconsistent  with  the  relief  depth  of  the  topography 
to  be  milled.  As  a  result,  milling  may  require  a  process  of  many  steps;  each  step  removing  a  thickness 
of  material  less  than  the  cutting-edge  length  of  the  tool.  A  process  based  on  use  of  a  single  tool  capable 
of  resolving  all  features  of  the  topography  would  be  too  inefficient  for  practical  use.  Instead,  a  multistep 
process  is  followed  in  which  the  size  of  the  tool  is  progressively  decreased  with  each  cutting  pass.  Such  a 
process  is  typically  categorized  into  two  stages  of  milling  operations:  roughing  and  finishing. 

In  the  roughing  stages,  one  or  more  flat  end  mills  are  used  to  remove  excess  material  from  areas  where 
no  detailed  profiling  of  the  surface  is  required.  Various  strategies  may  be  followed  in  preparing  the  surface 
for  finishing  [75].  One  class  of  method  removes  excess  material  by  machining  a  series  of  offset  surfaces 
of  increasing  fidelity.  Once  a  sufficient  amount  of  material  has  been  removed,  the  remaining  detailing  is 
completed  during  finishing  stages.  An  alternative  method,  used  here,  divides  the  surface  into  a  number  of 
equal- height  contours  and  removes  material  sequentially  from  each  layer.  The  resulting  terraced  surface  is 
sculpted  into  the  desired  form  during  the  finishing  stages. 

In  the  finishing  stages,  one  or  more  ball-nose  end  mills  are  used  to  cut  increasingly  detailed  surfaces. 
During  each  finishing  stage,  a  surface  derived  from  the  desired  topography,  but  having  a  level  of  detail 
appropriate  to  the  cutter  radius  used,  is  milled.  In  the  final  finishing  stage,  the  smallest  cutter,  having 
vj  <  p,  is  used  to  mill  the  desired  topography. 

In  order  to  machine  the  surface  accurately,  without  removing  excess  surface  material,  the  center  point 
of  the  spherical  tip  of  the  end  mill  must  be  offset  from  the  surface  by  w  in  the  direction  of  the  outward 
local  normal  vector  to  the  surface  [74].  Equivalently,  this  can  be  expressed  as  the  condition  that  the  center 
point  of  a  ball-nose  end  mill  of  tip  radius  w  moves  in  the  tu-evolute  surface  associated  with  the  desired 
rough  surface.  This  is  illustrated  in  Fig.  9.  Note  that  the  evolute  surface  allows  for  accurate  machining 
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only  where  w  <  p.  In  the  case  where  w  >  p  the  evolute  surface  is  no  longer  single  valued  and  gouging 
occurs.  Typically,  the  tool  path  consists  of  a  series  of  tool  passes:  closely  spaced  profiles  along  the  evolute 
surface,  which  the  cutter  is  made  to  raster  across  in  the  process  of  machining  the  surface.  The  final  fidelity 


FIG.  10:  A  schematic  representation  of  scalloping  caused  by  by  machining  a  surface  with  a  ball-nose  end  mill. 

of  the  finishing  stage  is  determined  not  only  by  the  smallest  tip-radius  of  ball-nose  end  mill  used,  but  also  by 
the  tool-pass  interval,  which  is  the  nominal  linear  spacing  between  adjacent  profiles  followed  by  the  cutter. 
An  ideal  ball- nose  end  mill  moving  on  an  ideal  evolute  surface  for  which  w  <  p  contacts  z( r)  only  on  the 
profiles  (tool  passes)  comprising  the  tool  path.  Consequently,  for  any  finite  tool  pass  interval,  the  finishing 
stage  leaves  excess  material  on  the  surface,  as  shown  in  Fig.  10.  The  height  of  these  so-called  scallops  is 
determined  by  the  tool-pass  interval,  the  tip  radius  of  the  ball-nose  end  mill,  and  the  local  topography  of 
z( r)  [73,  74].  For  a  flat  surface,  the  scallops  are  of  uniform  height  £  =  vo  —  \jw  —  (A/2)2,  where  A  is  the 
tool-pass  interval.  The  primary  requirement,  after  fidelity,  is  scallops  of  a  sufficiently  small  size  that  neither 
the  collective  effect  of  isolated  scattering  from  single  scallops  nor  the  effect  of  the  grating  scattering  resulting 
from  their  aggregate  are  significant  in  the  frequency  range  of  interest. 

An  example  of  a  rough  surface  physically  realized  using  the  techniques  described  here  is  shown  in  Figs. 
11-12.  In  Fig.  11,  the  desired  topography  is  shown  and  compared  with  the  topography  of  the  realized 
sample  as  measured  by  a  touch-trigger  probe,  shown  in  Fig.  12.  The  maximum  difference  between  these 
two  measured  topographies  is  <  1  mm. 
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FIG.  11:  Topography  of  a  computer-generated  realization  of  band-limited  FBM. 


C.  Verifying  Properties  of  Rough-Surface  Realizations 

Measurement  and  characterization  of  physical  realizations  of  band-limited  FBM  and  other  random  rough 
surfaces  are  needed  both  to  confirm  that  manufactured  physical  realizations  have  achieved  the  desired  form 
and  to  characterize  existing  rough  surfaces. 


1.  Measurement 

Reviews  of  techniques  for  measurement  of  rough-surface  topography  can  be  found  in  Refs.  76  and 
7.  Many  sophisticated  techniques  exist,  but  simple  contact  methods  are  sufficient  to  characterize  surface 
topography  over  length  scales  relevant  for  acoustical  scattering.  Such  methods  include  traditional  stylus 
methods  (profilometry)  as  described  in  Ref.  [76]  and  the  modern  and  conceptually  simple  kinematic-resistive 
touch-trigger  probe  used  in  coordinate  measurement  machines  (CMMs)  [77,  78,  79]. 

For  the  work  described  here,  a  touch-trigger  probe  was  used.  Unlike  profilometric  methods  that  provide 
a  continuous  reporting  of  surface  height,  touch-trigger  probes  specify  surface  coordinates  at  discrete  points. 
For  rough  surfaces  this  translates  to  specification  of  the  surface  height  at  discrete  points  on  a  uniform  two- 
dimensional  grid,  a  representation  termed  a  digital  elevation  model  (DEM).  This  representation  is  analogous 
to  that  used  for  the  numerical  generation  of  rough  surface  topography  described  in  Sec.  Ill  A,  and  therefore 
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FIG.  12:  Topography  of  a  physical  realization  manufactured  from  the  numerical  realization  shown  in  Fig.  11  as  measured  by 

a  touch-trigger  probe. 


has  the  same  large-/c  cutoff  described  in  Sec.  Ill  A  1. 

Similarly,  the  tip  radius  of  the  touch-trigger  probe  must  meet  the  criteria  given  in  Sec.  Ill  B  in  order  to 
accurately  represent  the  surface  topography.  If  the  tip  radius  w  does  not  satisfy  w  <  p,  the  touch-trigger 
probe  will  report  truncated  depths  for  those  valleys  where  w  >  p.  This  effect  of  the  spherical  radius  of  the 
probe  is  similar  to  that  of  a  low-pass  wavenumber  filter  applied  to  the  topography  [80,  81]. 

Data  obtained  by  physical  measurement  commonly  contain  errors  in  the  form  of  piston  or  tilt  (equivalent 
to  biases  of  c)  =  1  profiles  described  by  zeroth-  or  first-order  polynomials,  respectively)  due  to  misalignment 
between  the  surface  and  the  measuring  machine.  This  can  be  removed  by  fitting  a  plane  to  the  surface  and 
then  subtracting  this  plane  from  the  measured  DEM.  The  resulting  detrended  data  have  a  modified  PSD,  as 
described  in  Sec.  (II B).  Additionally,  detrending  is  appropriate  only  if  the  sample  is  large  enough  relative 
to  the  small- fc  cutoff  of  the  PSD  that  there  is  no  underlying  bias  in  the  topography  of  the  sample  itself. 


2.  Analysis 

Analysis  of  DEMs,  either  numerically  generated  or  measured  from  real  surfaces,  can  serve  a  number  of 
purposes.  For  numerically  derived  DEMs,  analysis  can  determine  whether  a  realization  is  representative  of 
the  ensemble  from  which  it  was  drawn.  For  physically  derived  DEMs,  analysis  can  characterize  and  unknown 
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physical  surface  or  validate  that  a  manufactured  surface  has  the  desired  properties. 

Rough-surface  topography  is  primarily  characterized  by  its  PSD.  Consequently,  analysis  of  DEMs  typ¬ 
ically  comprises  spectral  estimation  followed  by  parametric  descriptions  of  the  spectral  estimate,  such  as 
Wq  and  7t).  While  many  methods  of  spectral  estimation  exist  (see,  e.g.,  Ref.  82),  PSDs  of  power-law  form 
present  certain  complications  [83].  Spectral  estimates  based  on  the  Fourier  transform  of  a  finite-length  pro¬ 
file  suffer  from  spectral  leakage  due  to  convolution  in  k- space  of  the  Fourier  transform  of  the  topography 
with  the  Fourier  transform  of  the  windowing  function.  As  noted  by  Austin  et  al.  [83],  spectral  leakage  is  a 
particular  problem  for  PSDs  of  power-law  form  because  leakage  of  high-power  components  at  small  k  raises 
the  estimated  power  of  components  at  large  k  resulting  in  an  overall  bias  of  the  spectral  estimate  toward 
PSD  of  shallower  slope.  This  results  in  underestimation  of  7^  and  overestimation  of  w?,. 

Use  of  spatial  windows  having  rapidly  decaying  Fourier  transforms  can  partially  address  this,  but  only 
at  the  expense  of  decreased  spectral  resolution.  Alternatively,  DEMs  can  be  prewhitened  before  spectral 
estimation  in  order  to  reduce  spectral  leakage  and  the  estimated  PSD  corrected  by  a  multiplicative  factor 
to  account  for  prewhitening  [8,  13,  26,  83].  Prewhitening  is  easily  achieved  by  taking  the  first  differences  of 
the  data  [8,  26]. 

An  alternative  to  spectral  characterization  of  DEMs  is  characterization  in  terms  of  fractal  parameters  H 
and  r.  While  fractal  parameters  are  typically  estimated  from  PSD,  they  can  also  be  estimated  by  nonspectral 
methods.  This  is  necessary  in  the  case  of  DEM  because  their  inherent  band  limiting  causes  the  relationship 
between  H  and  7^  to  break  down  [49,  50,  53]. 

IV.  SUMMARY 

The  rough  topography  of  an  ocean  bottom  can  be  modeled  as  an  isotropic  two-dimensional  random  field 
having  a  power  spectral  density  of  power-law  form.  While  idealized  random  fields  having  power  spectral 
density  of  this  form  are  self-affine  stochastic  fractals  (realizations  of  FBM),  manufacturing  constraints  impose 
high-  and  low- A:  cutoffs  on  the  power-law  form  of  the  PSD.  As  a  result,  physically  realized  scale  models  of 
ocean  bottoms  are  approximately  self-affine  surfaces  (realizations  of  band-limited  FBM)  for  which  scaling 
properties  exist  only  over  a  finite  hierarchy  of  scales  and  conventional  relations  between  PSD  and  fractal 
geometry  break  down.  Rather  than  serve  as  a  complication,  the  breakdown  of  the  idealized  fractal  model 
produces  meaningful  relations  between  concepts  from  fractal  geometry  and  conventional  statistical  metrics. 
These  new  relations  explicitly  involve  wavenumber  cutoffs,  which,  as  shown  in  Sec.  II D,  are  as  influential 
on  the  appearance  of  fractal  surfaces  as  the  parameters  describing  the  PSD. 

Discretely  sampled  realizations  of  band-limited  FBM  can  be  numerically  generated  by  Fourier-transform 
methods,  with  high-  and  low- A:  cutoff  corresponding  to  sampling  interval  and  array  size,  respectively.  Physical 
realizations  of  these  numerical  representations  can  be  accurately  manufactured  by  CAM  techniques  using 
a  CNC  milling  machine  under  the  condition  that  the  tip  radius  of  the  ball- nose  end  mill  used  for  the  final 
finishing  is  smaller  than  the  minimum  positive  radius  of  curvature  of  the  surface. 

In  this  report  the  recent  concept  of  approximately  self-affine  profiles  and  surfaces  has  been  newly  intro¬ 
duced  into  the  context  of  modeling  the  rough  topography  of  the  ocean  bottom  and  the  fabrication  of  physical 
models.  Although  physical  realizations  of  FBM  have  been  previously  fabricated  [71,  72],  prior  work  has  not 
considered  the  effects  of  the  band-limiting  imposed  by  the  fabrication  process  on  the  fractal  properties  of 
the  resulting  physical  realization.  It  has  been  shown  that,  given  restrictions  in  bandwidth,  the  parameters 
of  the  power-law  spectrum  typically  used  in  the  field  of  underwater  acoustics  to  characterize  the  surface,  are 
no  longer  simply  related  to  the  more  fundamental  description  of  the  surface  in  terms  of  fractal  geometry. 
Moreover,  as  Sec.  II D  illustrates,  the  location  of  the  low -k  cutoff  can  markedly  affect  even  the  appearance 
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of  approximately  self- affine  band-limited  FBM.  Additionally,  this  report  has  given  algorithms  for  numerical 
generation  of  surface  topography  of  this  type  and  guidelines  for  CAM  that  are  specific  to  this  topography. 
The  discussion  concerning  numerical  generation  contributes  a  more  explicit  clarification  of  the  differences 
between  previously  published  spectral  methods. 
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APPENDIX  A:  REALIZATIONS  OF  BAND-LIMITED  FRACTAL  SURFACES 


This  appendix  contains  the  still  frames  of  the  embedded  movie  in  Fig.  3,  as  described  in  Table  I. 


FIG.  13:  Realization  of  a  band-limited  fractal  surface  with  102  =  1  •  10  5  m4,  72  =  2.0,  and  ka,L  =  0.05.  (See  Table  I.) 
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APPENDIX  B:  DETAILS  OF  ALGORITHM  FOR  GENERATION  OF  SURFACE 
REALIZATIONS 

Because  z( r)  is  real,  Z( k)  must  be  conjugate  symmetric.  Efficient  computation  of  z( r)  from  Z( k)  by  the 
DFT  using  the  FFT  algorithm  requires  that  the  matrices  have  dimensions  that  are  integer  powers  of  two. 
A  conjugate  symmetric  M  by  N  matrix  where  both  M  and  N  are  even  is  of  the  form 
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where  the  subscripts  indicate  the  row  and  column,  respectively.  In  those  cases  where  the  indices  appear  to  be 
incorrect,  it  simply  represents  the  conjugate  symmetry  of  the  matrix.  The  matrix  of  uniformly  distributed 


random  numbers  a  used  in  Eqs.  (42  and  44)  must  have  odd  symmetry,  as  illustrated  by  the  antisymmetric 
matrix  below 

al,l  al,N  al,N-l  ' 

ai,f+2  “l.f+l 

’  *  ’  al,AT 

al,2  ~aM,N  aM,N  —  l  * 

••  “aM,f+2  1 

af+i,i  : 

1  “f  +  l.f+l 

...  :  .  (B2) 

—  CL  M 

2 

:  +  l  af+2,f+2 

~a2,l  a2,N  ~a2,N-l  ' 

_ a2,f+2  _ a2,f+l  aM,f+2 

'  '  '  aM,N 

In  contrast,  the  matrix  of  uniformly  distributed  random  numbers  f3  used  in 
symmetry,  as  illustrated  by  the  symmetric  matrix  below 

Eq.  (44)  must  have  even 

al,l  al,N  al,N—l  '  ' 

'  °l,f+2  ai,f+l  . 

'  al,N 

al,2  aM,N  aM,N-l  *  * 

'  °M,f+2  : 

af+i,i  : 

:  “f+l.f  +  l  . 

■  :  .  (B3) 

CL  m 

2  > 

:  “f.f  +  l  ®f+2,f+2  " 

a2,l  a2,N  a2,N-l  ’  * 

°2,f^+2  a2,f^+l  aM,^f+2 

'  aM,N 
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APPENDIX  C:  MATLAB  FUNCTIONS 


This  appendix  contains  the  header  portions  of  matlab  functions  used  to  numerically  generate  and 
characterize  band-limited  FBM.  The  full  functions  along  with  required  support  functions  are  included  on 
this  CD-ROM.  The  matlab  functions  were  developed  under  matlab  6.1  R12  and,  in  some  cases,  require 
the  Signal  Processing  Toolbox  (or  equivalent  functions  written  by  the  user)  for  operation.  All  software  on 
this  CD-ROM  is  distributed  freely,  but  is  subject  to  the  copyrights  and  terms  of  use  specified  in  the  headers 
of  the  individual  files. 


1.  Numerical  generation  of  band-limited  FBM  surfaces 


function  surf _z=generate_surf ace (varargin) 

l  □  =generate_surf ace (f rq, c , Lx , Ly , dx , dy , spectratype , w_2_gamma_2 , 

%  LF_param ,pname ,randstate) 

/ 

1  Generates  a  2-D  random  rough  surface  with  a  power-law  spectrum 
“/,  and  low-frequency  roll-off.  The  surface  and  parameters  are  saved 
"/.  to  a  file.  Plots  of  the  surface  can  be  displayed. 

°/o 

°/0  If  called  as  a  function  with  input  parameters,  the  parameters 

"/.  specified  in  the  function  call  are  used  instead  of  the 

“/,  parameters  specified  in  the  file  and  plotting  is  disabled. 

1 

1  J.E.  Summers,  2005 


%  Reference 

mmmmmmmmmmmmmmmmmmxmmmxm 

°/= 


"/.  J.  E.  Summers,  R.  J.  Soukup,  and  R.  F.  Gragg,  "Characterization 
“/,  and  fabrication  of  synthetic  rough  surfaces  for  acoustical  scale 
“/,  -model  experiments,"  NRL/MR-MM/7140 — 05-8871,  U.  S.  Naval  Research 
*/.  Laboratory,  Washington,  D.C.,  2005 
% 


rararararararaffirararanmnraffimmraffiraffix 


mmmmmmmmmmmmmmmmmxmmmxmm 


mmmmmMmmmmmmmmmmmmxmiimmi! 

1  Dependencies 


%  requires : 

*/.  Signal  Processing  Toolbox:  fftshift.m 

“/„  calcspectra.m  :  computes  the  2-D  roughness  spectrum 

*/.  rand2.m  :  computes  the  2-D  realization 

“/.  odd.m  :  rounds  to  the  nearest  odd  number  (time-frequency  toolbox) 


1  can  be  used  with: 

%  characterize_surf ace .m  :  characterizes  the  surface  and  generates 
“/.  a  text  file  of  descriptors 

%  export_surface .m  :  exports  the  surface  geometry  into  common 
*/.  CAD/CAM  formats  (*.asc,  *.stl) 


ffiffirararamfflrarafflnmnraffimmraramrara 


mmmiimmmmmmx  nmmnmmmmmxmmmmii 

%  Input  Parameters 


if  nargin<l 

*/.  target  measurement  frequency  (Hz) 
f rq  =  300000; 


l  Conditions  and  terms  for  the  use  of  this  software 

1 


°/0  This  software  may  be  used,  reproduced,  modified,  and  distributed 
l  by  the  user  for  noncommercial  purposes  but  is  Copyright  (C) 
l  Jason  E.  Summers  (summers@abyss.nrl.navy.mil). 

°/= 


7.  Users  have  the  right  to  use  the  software  with  the  following  terms 
“/,  and  conditions: 


*/. 


%  Users  agree  that  any  copies  of  the  software  will  contain  the  same 
%  proprietary  notices  and  warranty  disclaimers  which  appear  in  this 
%  software . 


% 

%  Users  can  make  modifications  to  this  software  but  must 
%  include  references  to  the  original  software  and  author. 

% 

%  If  a  user  makes  any  modifications  to  this  software  and  the 
%  modified  software  is  redistributed,  then  the  user  must  notify  the 
%  author  of  the  original  software  by  an  e-mail  sent  to 
%  summers@abyss.nrl.navy.mil. 

1 


*/.  The  original  author  shall  be  credited  should  this  software  be  used 
%  in  any  form  or  written  about  in  any  publication. 

°/= 


*/.  This  software  is  provided  "as  is",  without  warranty  by  the  author 
°/0  or  by  the  United  States  Naval  Research  Laboratory  (NRL)  .  Although 
*/.  the  software  has  been  tested  on  Windows  machines  under 
%  MATLAB  6.1  R12,  neither  the  author  nor  NRL  make  any  warranties 

%  relating  to  the  software’s  performance  on  this  or  any  other 

%  platforms.  The  author  and  NRL  make  no  express,  implied,  or 

%  statuatory  warranty  of  any  kind  for  this  software  including,  but 

%  not  limited  to,  any  warranty  of  performance,  merchantability,  or 

%  fitness  for  a  particular  purpose.  In  no  event  shall  the  author  or 
%  NRL  be  liable  to  the  user  or  any  thrid  party  for  any  loss  or  for 
•/.  any  indirect,  special,  punitive,  examplary,  incidental,  or 
%  consequential  damages  arising  from  the  use,  possession  or 
°/0  performance  of  this  software,  even  if  the  author  or  NRL  have 
•/.  been  advised  of  the  possibility  thereof. 

°/o 


%  Any  comments  on  errors,  or  hints  about  possible  improvements  or 
•/.  extensions  are  appreciated. 

mmmmmmmmmnmmmmmmmmnmmmm 


%  speed  of  sound  (meters/sec) 

*/.  GPM  RSS2 :  1482  m/s 
c  =  1482; 

%  desired  surface  dimensions  (meters) 

Lx  =  1; 

Ly  =  l; 

%  desired  spacing  between  sample  points  (meters) 
dx  =  le-3; 
dy  =  le-3; 


%  form  of  the  spectrum 
spectratype=’alg’ ; 

*/.  *  ’moe’:  Moe’s  spectrum  (Gaussian  LF  cutoff)  [Eq.  (15)] 

*/.  *  ’pwr’  :  power  law  [Eq.  (10)] 

%  *  ’alg’ :  power  law  with  algebraic  LF  cutoff  [Eq.  (13)] 

y.  *  ’box’  :  power  law  with  box-car  HF  and  LF  cutoff  [Eq.  (12)] 

*/.  *  ’dif’:  power  law  with  LF  cutoff  expressed  as  the 

%  difference  of  two  pure  power  laws  [Eqs.  (16)] 

%  for  the  realization  either  the  inverse  transform  of  the  roughness 
%  spectrum  can  be  done  by  assuming  random  phases  with  a  Gaussian 
*/.  distribution  only  (amprand=0)  or  by  adding  a  random  amplitude 
%  factor  given  by  a  Rayleigh  distribution  (amprand=l) 
amprand=l ; 

%  spectral  strength  (meters ~4) 
w_2  =  2 . Oe-5 ; 

%  two-dimensional  spectral  exponent 
%  (dimensionless  number  between  2  and  4) 
gamma_2  =  3.5; 

%  low-frequency  cutoff  parameter 
LF_param  =  20; 

*/.  try  0.05  for  a  in  moe  (knee  at  K~30;  scales  in  K  linearly  with 
%  l/LF_param) 

*/.  try  20  for  mod 

“/path  name  in  which  to  save  created  surface 
pname= ’ D : \NRL\ ’ ; 

%  number  of  the  random  iteration  to  create 
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*/.  (0  to  2~1492,  which  is  the  period  of  Matlab’s  internal 

*/.  random  number  generator) 

randstate=333; 

7  turn  on  and  off  plotting 

plotsurf=0;  7  1  to  show  plots,  0  to  not 

plotspec=0;  7  1  to  show  plots  of  spectrum,  0  to  not 

elseif  nargin==12 


w_2=varargin{8} ; 
gamma_2=varargin{9} ; 

LF_param=varargin{10} ; 
pname=varargin{ll} ; 
randstate=varargin{12} ; 

else 

error ( ’generate_surf ace :  wrong  number  of  input  arguments’) 


f rq=varargin{l} ; 
c=varargin{2} ; 
Lx=varargin{3} ; 
Ly=varargin{4} ; 
dx=varargin{5} ; 
dy=varargin{6} ; 
spectrtype=varargin{7} ; 


end 


ffiffirararamfflrarafflnmnrammmraramrara 


2.  Characterization  of  discrete-sampled  rough  surfaces 


function  varargout=characterize_surf ace (varargin) 

7.  [w_2 , gamma_2]  =characterize_surf ace  (Spname ,  Sf name) 

7. 

7.  Where  Spname  and  Sfname  are  the  path  name  and  file  name  of  a  file 
7.  created  by  generate_surf ace .m. 

7 

7  or 
7. 

7o  [w_2 , gamma_2]  =characterize_surf ace  (c , dx , dy ,  surf _z , f rq) 

7. 

7.  Where  c  is  the  sound  speed,  dx  is  sample  spacing  in  x,  dy  is  the 
7.  sample  spacing  in  y,  and  surf_z  is  a  matrix  of  surface  heights. 

7. 

7.  or 

7. 

7.  [w_2 , gamma_2] =characterize_surf ace (c , surf _x , surf _y , surf _z , f rq) 

7. 

t  Where  c  is  the  sound  speed,  surf_x  is  a  matrix  of  x  coordinates, 

*/.  surf_y  is  a  matrix  of  y  coordinates,  and  surf_z  is  a  matrix  of 
7.  surface  heights. 

7. 

7.  Calculates  parameters  of  a  2D  random  rough  surface  with  a  power-law 
7.  spectrum,  which  are  then  displayed  on  the  screen  and  output  to  a 
*/.  text  file.  The  following  parameters  are  calculated: 

7. 

/  max  and  min  surface  height 
7.  rms  surface  height 
7.  minimum  curvature 

7.  smallest  ball-nose  end  mill  required  to  mill  surface  accurately 
7.  correlation  length 
7.  w_2 
7.  gamma_2 

7.  fractal  dimension 
7.  Hurst  exponent 
7.  topothesy 

7. 

7.  If  called  as  a  function,  plotting  and  screen  displays  are  disabled. 

7. 

7.  If  called  with  no  outputs,  results  are  displayed  on  screen  and 
7«  plotted. 

7. 

7.  If  called  with  no  inputs,  a  default  file  is  opended  and  analyzed. 

7. 
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7.7.7.7.7.7.7.7X/.7.7.7.7.7.7X/.7.7.7.7.7.7X/.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7. 

7.  Reference 

7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.77.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7 

7. 


7.  J.  E.  Summers,  R.  J.  Soukup,  and  R.  F.  Gragg,  "Characterization 
7.  and  fabrication  of  synthetic  rough  surfaces  for  acoustical  scale 
7.  -model  experiments,"  NRL/MR-MM/7140 — 05-8871,  U.  S.  Naval  Research 
7.  Laboratory,  Washington,  D.C.,  2005 

7. 


7.7.7.77.7.7.7.7.7.7.7.77.7.7.7.7.7.7.7.77.7.7.7.7.7.7.7.77.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.7.77.7.7.7 


7.7.7.7777.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.7777.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.7.777.7.7.7. 

7.7.7.7777.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.7777.7.7.7.7.7.777.7.7.7.7.7.777.7.7: 

7.  Conditions  and  terms  for  the  use  of  this  software 

7. 


7.  This  software  may  be  used,  reproduced,  modified,  and  distributed 
7.  by  the  user  for  noncommercial  purposes  but  is  Copyright  (C) 

7.  Jason  E.  Summers  (summers@abyss.nrl.navy.mil). 

7. 


7.  Users  have  the  right  to  use  the  software  with  the  following  terms 
7.  and  conditions: 


7. 


7  Users  agree  that  any  copies  of  the  software  will  contain  the  same 
7  proprietary  notices  and  warranty  disclaimers  which  appear  in  this 
7  software . 

7 

7  Users  can  make  modifications  to  this  software  but  must 
7  include  references  to  the  original  software  and  author. 

7 

7  If  a  user  makes  any  modifications  to  this  software  and  the 
7  modified  software  is  redistributed,  then  the  user  must  notify  the 
7  author  of  the  original  software  by  an  e-mail  sent  to 
7  summers@abyss.nrl.navy.mil. 

7 

7  The  original  author  shall  be  credited  should  this  software  be  used 
7  in  any  form  or  written  about  in  any  publication. 

7 

7.  This  software  is  provided  "as  is",  without  warranty  by  the  author 
7  or  by  the  United  States  Naval  Research  Laboratory  (NRL) .  Although 
7  the  software  has  been  tested  on  Windows  machines  under 
7  MATLAB  6.1  R12,  neither  the  author  nor  NRL  make  any  warranties 

7  relating  to  the  software’s  performance  on  this  or  any  other 

7  platforms.  The  author  and  NRL  make  no  express,  implied,  or 

7  statuatory  warranty  of  any  kind  for  this  software  including,  but 

7  not  limited  to,  any  warranty  of  performance,  merchantability,  or 

7  fitness  for  a  particular  purpose.  In  no  event  shall  the  author  or 

7  NRL  be  liable  to  the  user  or  any  thrid  party  for  any  loss  or  for 

7  any  indirect,  special,  punitive,  examplary,  incidental,  or 
7  consequential  damages  arising  from  the  use,  possession  or 
7  performance  of  this  software,  even  if  the  author  or  NRL  have 
7  been  advised  of  the  possibility  thereof. 

7 


7.  Any  comments  on  errors,  or  hints  about  possible  improvements  or 
7.  extensions  are  appreciated. 

7.7.7.7.7.7.7.7.7.7.77.7.7.7.7.7.7.7.7.7.7.7.7.7.777.7.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.7.77.7.7.7.7.7.7.777.7.7.7.7.7.7 


7.7.7.7.7.77777.7.7.7.77.777.7.7.7.77777.7.7.7.777777.7.7.7.7.7.7.77.7.7.7.7.7.777.7.7.7.777777.7.7.7.77777. 


7.7.7.7.7.77.7.7.7.7.7.7.77777.7.7.7.7.7.777.7.7.7.77.777.7.7.7.7.77.7.7.7.7.7.7.7.7.777.7.7.7.777777.7.7.7.7.7777. 
7.  Dependencies 

777.7.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7 


7  requires : 

7  Signal  Processing  Toolbox:  xcorr.m 
7  linfit.m  :  a  least-squares  linear-fitting  routine 
7  calc_curve.m  :  calculates  curvature  of  a  surface 

7 


7  can  be  used  with: 

7  generate_surface  .m 


7  export_surf ace .] 


7 


:  generates  a  2D  random  rough  surface 
with  a  specified  power-law  spectrum 
exports  the  surface  geometry  into  common 
CAD/CAM  formats  (*.asc,  *.stl) 


777.7.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.7.777.7.7.7.7.7.7.7.7.7.7.7.7.7.777.7.7.7.7.7.7 


if  nargout==l  |  nargout>2 

error (’ characterize_surf ace :  wrong  numer  of  output  arguments’) 


end 


77.7.7.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.7.77.7.7.7.7.7.7.7.77.7.7.7.7.7.77.7.7.7.7.7.7.777 

7  Input  Parameters 

777.7.7.7.7.7.7.7.77.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7.7.777.7.7.7.7 


if  nargin==0 

7.  file  name  containing  the  surface  realization  to  be  characterized 
Spname=’D: \NRL\GeoacousticPhysicalModeling\2004\ ’ ; 

Sf name= ’ W2_le-005_G2_3_rP_333 ’ ; 

7.  load  the  surface  realization  data 
load ([Spname  Sfname  ’ .mat’] , ’-mat’) 

7.  control  plotting 

plotstats=l;  7.  1  to  show  plots  of  statistics,  0  to  not 
plotf its=0; 
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normalize_plot=0 ;  %  1  to  normalize  the  S  vs .  k  plots,  0  to  not 

elseif  nargin==2 

if  ischar (varargin{l}) 

Spname=varargin{l} ; 

Sf name=varargin{2} ; 

1  load  the  surface  realization  data 
load( [Spname  Sfname  ’ .mat’] , ’-mat’) 

end 


dx=varargin{2} ; 
dy=varargin{3} ; 

else 

surf _x=varargin{2} ; 
surf _y=varargin{3} ; 

end 

end 

else 

error ( ’ characterize_surf ace:  wrong  number  of  input  arguments’) 


elseif  nargin==4 


end 


if  isnumeric(varargin{l}) 
c=varargin{l}; 
surf _z=varargin{4} ; 
if  size(varargin{2})== [1  1] 
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GLOSSARY  OF  SYMBOLS 


CAM 

Computer-aided  machining. 

CMM 

Coordinate  measurement  machine. 

DEM 

Digital  elevation  model. 

PDF 

Power  density  function. 

PSD 

Power  spectral  density. 

<>e 

Ensemble  average, 

<>S 

Spatial  average. 

71 

One-dimensional  roughness  exponent. 

72 

Radial  roughness  exponent. 

A 

Tool-pass  interval. 

Kr 

Mean  curvature. 

Kk 

Gaussian  curvature. 

«+ 

Maximum  principal  curvature. 

K— 

Minimum  principal  curvature. 

t 

Scallop  height. 

w 

Radius  of  ball-nose  end  mill. 

p 

Minimum  positive  radius  of  curvature 

p(r) 

Radius  of  curvature  of  z( r)  at  r. 

a 

Root-mean-square  roughness. 

T 

Topothesy. 

B 

Box  dimension. 

C 

Autocorrelation. 

Co 

Auto  covariance. 

d 

Euclidean  dimension. 

d 

Topological  dimension. 

fz(0 

Probability  density  function  of  z. 
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ho  Reference  length. 

H  Hausdorff-Besicovitch  dimension. 

Jo  Zeroth-order  Bessel  function. 

k  Radial  spatial  wave  number. 

Upper  cutoff  of  k. 
k _  Lower  cutoff  of  k. 

kn  Upper  cutoff  of  k. 

kL  Lower  cutoff  of  k. 

kp  Location  of  spectral  peak. 

kx  One- dimensional  spatial  wave  number. 

kx  Minimum  resolvable  kx. 

kx  Maximum  resolvable  kx. 

kciL  Lower  cutoff  of  k. 

i  Correlation  length. 

P(k)  Radial  power  spectral  density. 

p(kx)  Linear  power  spectral  density. 

r  Two-dimensional  vector  of  position  with  coordinates  (x,y). 

R  Translation  |ri  —  7*2 1. 

W1  Set  of  real  number  in  n  dimensions. 

S( R)  Structure  function. 

wi  One-dimensional  spectral  strength. 

W2  Radial  spectral  strength. 

z( r)  Surface  topography  (height)  as  a  function  of  position. 

z0(r)  Deterministic  large-scale  roughness. 

Z( k)  Wavenumber  (k- space)  representation  of  the  surface. 

<  z  >  Mean  surface  height. 
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